ESPResSo
Extensible Simulation Package for Research on Soft Matter Systems
Loading...
Searching...
No Matches
EKReaction.hpp
Go to the documentation of this file.
1/*
2 * Copyright (C) 2022-2023 The ESPResSo project
3 *
4 * This file is part of ESPResSo.
5 *
6 * ESPResSo is free software: you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation, either version 3 of the License, or
9 * (at your option) any later version.
10 *
11 * ESPResSo is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program. If not, see <http://www.gnu.org/licenses/>.
18 */
19
20#pragma once
21
22#include "config/config.hpp"
23
24#ifdef WALBERLA
25
26#include "EKReactant.hpp"
27#include "LatticeIndices.hpp"
28#include "LatticeWalberla.hpp"
29
33
36
38
39#include <algorithm>
40#include <memory>
41#include <sstream>
42#include <stdexcept>
43#include <string>
44#include <vector>
45
47
48class EKReaction : public AutoParameters<EKReaction, LatticeIndices> {
49public:
50 [[nodiscard]] std::shared_ptr<::walberla::EKReactionBase>
51 get_instance() const {
52 return m_ekreaction;
53 }
54
58
59protected:
60 auto get_agrid(VariantMap const &args) const {
61 auto lattice = get_value<std::shared_ptr<LatticeWalberla>>(args, "lattice");
62 return get_value<double>(lattice->get_parameter("agrid"));
63 }
64
66 auto const tau = get_value<double>(args, "tau");
67 auto const agrid = get_agrid(args);
68 auto reactant = get_value<std::vector<Variant>>(args, "reactants");
69
70 auto get_order = [](Variant const &v) {
71 return get_value<double>(
72 get_value<std::shared_ptr<EKReactant>>(v)->get_parameter("order"));
73 };
74 auto const sum_alphas =
75 std::accumulate(reactant.begin(), reactant.end(), 0.,
76 [get_order](double sum, auto &element) {
77 return sum + get_order(element);
78 });
79
80 return tau / std::pow(Utils::int_pow<3>(agrid), sum_alphas - 1.);
81 }
82
83 template <typename F>
84 auto make_instance(VariantMap const &args, F &allocator) const {
85 auto lattice = get_value<std::shared_ptr<LatticeWalberla>>(args, "lattice");
86 auto reactants = get_value<std::vector<Variant>>(args, "reactants");
88 auto get_instance = [](Variant const &v) {
90 };
91 std::transform(reactants.begin(), reactants.end(), output.begin(),
93
94 auto const coefficient =
96
97 return allocator(lattice->lattice(), output, coefficient);
98 }
99
100 std::shared_ptr<::walberla::EKReactionBase> m_ekreaction;
102};
103
105public:
107 add_parameters({{"coefficient",
108 [this](Variant const &v) {
109 get_instance()->set_coefficient(
111 },
112 [this]() {
113 return get_instance()->get_coefficient() /
115 }}});
116 }
117
122};
123
125public:
128 {{"coefficient",
129 [this](Variant const &v) {
130 get_instance()->set_coefficient(get_value<double>(v) *
132 },
133 [this]() {
134 return get_instance()->get_coefficient() /
136 }},
137 {"shape", AutoParameter::read_only, [this]() {
138 return get_instance()->get_lattice()->get_grid_dimensions();
139 }}});
140 }
141
142 void do_construct(VariantMap const &args) override {
143 auto const agrid = get_agrid(args);
145 m_ekreaction_impl =
147 m_ekreaction = m_ekreaction_impl;
148 }
149
150 [[nodiscard]] Variant do_call_method(std::string const &method,
151 VariantMap const &parameters) override {
152 if (method == "set_node_is_boundary") {
153 auto const index = get_mapped_index(
155 get_instance()->get_lattice()->get_grid_dimensions());
156 m_ekreaction_impl->set_node_is_boundary(
157 index, get_value<bool>(parameters, "is_boundary"));
158 return none;
159 }
160 if (method == "get_node_is_boundary") {
161 auto const index = get_mapped_index(
163 get_instance()->get_lattice()->get_grid_dimensions());
164 auto const result = m_ekreaction_impl->get_node_is_boundary(index);
165 return Utils::Mpi::reduce_optional(context()->get_comm(), result);
166 }
167 return {};
168 }
169
170private:
171 std::shared_ptr<::walberla::EKReactionBaseIndexed> m_ekreaction_impl;
172};
173
174} // namespace ScriptInterface::walberla
175
176#endif // WALBERLA
Bind parameters in the script interface.
Variant get_parameter(const std::string &name) const final
void add_parameters(std::vector< AutoParameter > &&params)
Utils::Vector3i get_mapped_index(Utils::Vector3i const &index, Utils::Vector3i const &shape) const
Context * context() const
Responsible context.
void do_construct(VariantMap const &args) override
void do_construct(VariantMap const &args) override
Variant do_call_method(std::string const &method, VariantMap const &parameters) override
auto get_conversion_coefficient() const noexcept
std::shared_ptr<::walberla::EKReactionBase > m_ekreaction
auto get_agrid(VariantMap const &args) const
std::shared_ptr<::walberla::EKReactionBase > get_instance() const
auto calculate_bulk_conversion_factor(VariantMap const &args) const
auto make_instance(VariantMap const &args, F &allocator) const
std::vector< std::shared_ptr< EKReactant > > reactants_type
This file contains the defaults for ESPResSo.
T get_value(Variant const &v)
Extract value of specific type T from a Variant.
std::unordered_map< std::string, Variant > VariantMap
Definition Variant.hpp:69
boost::make_recursive_variant< None, bool, int, std::size_t, double, std::string, ObjectRef, Utils::Vector3b, Utils::Vector3i, Utils::Vector2d, Utils::Vector3d, Utils::Vector4d, std::vector< int >, std::vector< double >, std::vector< boost::recursive_variant_ >, std::unordered_map< int, boost::recursive_variant_ >, std::unordered_map< std::string, boost::recursive_variant_ > >::type Variant
Possible types for parameters.
Definition Variant.hpp:67
constexpr const None none
None-"literal".
Definition Variant.hpp:50
T reduce_optional(boost::mpi::communicator const &comm, std::optional< T > const &result)
Reduce an optional on the head node.
std::shared_ptr< EKReactionBase > new_ek_reaction_bulk(std::shared_ptr< LatticeWalberla > const &lattice, typename EKReactionBase::reactants_type const &reactants, double coefficient)
std::shared_ptr< EKReactionBaseIndexed > new_ek_reaction_indexed(std::shared_ptr< LatticeWalberla > const &lattice, typename EKReactionBase::reactants_type const &reactants, double coefficient)
static constexpr const ReadOnly read_only