8#include <compile_time_options.h>
12#include <deal.II/base/mpi.h>
13#include <deal.II/base/parameter_acceptor.h>
15#include <boost/signals2.hpp>
24 inline const std::string
dave =
25 "\nDave, this conversation can serve no purpose anymore. Goodbye.\n\n";
35 : ParameterAcceptor(
"B - Equation")
38 add_parameter(
"dimension", dimension_,
"The spatial dimension");
39 add_parameter(
"equation", equation_,
"The PDE system");
41 time_loop_executed_ =
false;
52 dave +
"No equation has been registered. Consequently, "
53 "there is nothing for us to do.\n"));
62 template <
typename Callable>
75 void dispatch(
const std::string ¶meter_file,
const MPI_Comm &mpi_comm)
77 ParameterAcceptor::prm.parse_input(parameter_file,
82 AssertThrow(dimension_ >= 1 && dimension_ <= 3,
83 dealii::ExcMessage(
dave +
84 "The dimension parameter needs to be "
85 "either 1, 2, or 3, but we encountered »" +
86 std::to_string(dimension_) +
"«\n"));
90 dave +
"No equation has been registered. Consequently, "
91 "there is nothing for us to do.\n"));
94 dimension_, equation_, parameter_file, mpi_comm, time_loop_executed_);
96 AssertThrow(time_loop_executed_ ==
true,
97 dealii::ExcMessage(
dave +
98 "No equation was dispatched "
99 "with the chosen equation parameter »" +
100 equation_ +
"«.\n"));
107 template <
typename Callable>
130 boost::signals2::signal<void(
int ,
131 const std::string & ,
132 const std::string & ,
152 std::string equation_;
156 bool time_loop_executed_;
165 template <
typename Description,
int dim,
typename Number>
167 bool write_detailed_description)
178 auto &prm = dealii::ParameterAcceptor::prm;
179 prm.enter_subsection(
"B - Equation");
180 prm.set(
"dimension", std::to_string(dim));
181 prm.set(
"equation", name);
182 prm.leave_subsection();
184 std::string base_name = name;
185 std::replace(base_name.begin(), base_name.end(),
' ',
'_');
186 base_name +=
"-" + std::to_string(dim) +
"d";
188 if (dealii::Utilities::MPI::this_mpi_process(MPI_COMM_SELF) == 0) {
189 const auto full_name =
190 "default_parameters-" + base_name +
"-description.prm";
191 if (write_detailed_description)
192 prm.print_parameters(full_name,
193 dealii::ParameterHandler::OutputStyle::PRM);
195 const auto short_name =
"default_parameters-" + base_name +
".prm";
196 prm.print_parameters(short_name,
197 dealii::ParameterHandler::OutputStyle::Short);
203 dealii::ParameterAcceptor::clear();
211 template <
typename Description,
typename Number>
216 std::cout <<
"Dispatch<Description, Number>::Dispatch() for »" << name
221 create_prm_files<Description, 1, Number>(name,
false);
222 create_prm_files<Description, 2, Number>(name,
true);
223 create_prm_files<Description, 3, Number>(name,
false);
227 [name](
const int dimension,
228 const std::string &equation,
229 const std::string ¶meter_file,
230 const MPI_Comm &mpi_comm,
231 bool &time_loop_executed) {
232 if (equation != name)
235 if (dealii::Utilities::MPI::this_mpi_process(mpi_comm) == 0) {
236 std::cout <<
"[INFO] dispatching to driver »" << equation
237 <<
"« with dim=" << dimension << std::endl;
240 AssertThrow(time_loop_executed ==
false,
243 "Trying to execute more than one TimeLoop object "
244 "with the given equation parameter »" +
247 if (dimension == 1) {
249 dealii::ParameterAcceptor::initialize(parameter_file);
251 time_loop_executed =
true;
252 }
else if (dimension == 2) {
254 dealii::ParameterAcceptor::initialize(parameter_file);
256 time_loop_executed =
true;
257 }
else if (dimension == 3) {
259 dealii::ParameterAcceptor::initialize(parameter_file);
261 time_loop_executed =
true;
void dispatch(const std::string ¶meter_file, const MPI_Comm &mpi_comm)
static void register_create_parameter_files(const Callable &callable)
static void create_parameter_files()
static void register_dispatch(const Callable &callable)
void create_prm_files(const std::string &name, bool write_detailed_description)
Dispatch(const std::string &name)
boost::signals2::signal< void(int, const std::string &, const std::string &, const MPI_Comm &, bool &)> dispatch
boost::signals2::signal< void()> create_parameter_files