Dune::Copasi 2.1.0
make_step_operator.hh
Go to the documentation of this file.
1#ifndef DUNE_COPASI_MODEL_MAKE_DIFFUSION_REACTION_STEP_OPERATOR_HH
2#define DUNE_COPASI_MODEL_MAKE_DIFFUSION_REACTION_STEP_OPERATOR_HH
3
6
7#if HAVE_METIS && DUNE_COPASI_CONCURRENT_ASSEMBLY
8#include <dune/pdelab/common/partition/metis.hh>
9#endif
10
11#include <dune/pdelab/common/partition/simple.hh>
12#include <dune/pdelab/common/execution.hh>
13#include <dune/pdelab/operator/operator.hh>
14
15#include <spdlog/spdlog.h>
16
17#include <memory>
18
20
21template<class LocalBasisTraits,
22 class Coefficients,
23 class Residual,
24 class ResidualQuantity,
25 class TimeQuantity,
26 PDELab::Concept::Basis Basis,
27 Dune::Concept::GridView CellDataGridView = typename Basis::EntitySet,
28 class CellDataType = double>
29[[nodiscard]] inline static std::unique_ptr<PDELab::Operator<Coefficients, Coefficients>>
30make_step_operator(const ParameterTree& config,
31 const Basis& basis,
32 std::size_t halo,
33 std::shared_ptr<const FunctorFactory<Basis::EntitySet::dimension>> functor_factory,
34 std::shared_ptr<const CellData<CellDataGridView, CellDataType>> grid_cell_data = nullptr)
35{
36
37 std::unique_ptr<PDELab::Operator<Coefficients, Coefficients>> one_step;
38 const auto& assembly_cfg = config.sub("assembly");
39
40 auto make_one_step_op = [&, functor_factory, grid_cell_data]<class ExecutionPolicy, PDELab::Concept::Basis OperatorBasis>(
41 ExecutionPolicy execution_policy, const OperatorBasis& operator_basis) {
42 spdlog::info("Creating mass/stiffness local operator");
43 const auto& time_step_cfg = config.sub("time_step_operator");
44 using LocalOperator =
45 LocalOperator<OperatorBasis, LocalBasisTraits, CellDataGridView, CellDataType, ExecutionPolicy>;
46 LocalOperator const stiff_lop(operator_basis,
48 config,
49 functor_factory,
50 grid_cell_data,
51 execution_policy);
52 LocalOperator const mass_lop(operator_basis,
54 config,
55 functor_factory,
56 grid_cell_data,
57 execution_policy);
58 return Dune::Copasi::make_step_operator<Coefficients, Residual, ResidualQuantity, TimeQuantity>(
59 time_step_cfg, operator_basis, mass_lop, stiff_lop);
60 };
61
62 auto entity_set = basis.entitySet();
63 std::size_t part_patches =
64 assembly_cfg.get("partition.patches", static_cast<std::size_t>(entity_set.size(0) / 40));
65
66#if DUNE_COPASI_CONCURRENT_ASSEMBLY
67 auto make_concurrent_one_step_op =
68 [&]<class ExecutionPolicy>(const ExecutionPolicy execution_policy) {
69 auto entity_set = basis.entitySet();
70 auto default_part_type =
71#if HAVE_METIS
72 "METIS";
73#else
74 "Simple";
75#endif // HAVE_METIS
76 auto part_type = assembly_cfg.get("partition.type", default_part_type);
77 auto part_coloring = assembly_cfg.get("partition.coloring", "None");
78 spdlog::info(" Concurrent Partitioning:");
79 spdlog::info(" Type: {}", part_type);
80 spdlog::info(" Patches: {}", part_patches);
81 spdlog::info(" Halo: {}", halo);
82 spdlog::info(" Coloring: {}", part_coloring);
83 if (part_type == "METIS") {
84#if HAVE_METIS
85 if (part_coloring == "None") {
86 return make_one_step_op(
87 execution_policy,
88 basis.subSpace(
89 PDELab::EntitySetPartitioner::Metis{ entity_set, part_patches, halo },
90 TypeTree::treePath()));
91 } else if (part_coloring == "DSatur") {
92 return make_one_step_op(
93 execution_policy,
94 basis.subSpace(PDELab::EntitySetPartitioner::MetisColored{ entity_set, part_patches, halo },
95 TypeTree::treePath()));
96 } else {
97 throw format_exception(IOError{}, "Not known coloring algorithm '{}' known", part_coloring);
98 }
99#else
100 throw format_exception(IOError{}, "'METIS' partitioner is not available on this executable");
101#endif // HAVE_METIS
102 } else if (part_type == "Simple") {
103 if (part_coloring == "None") {
104 return make_one_step_op(
105 execution_policy,
106 basis.subSpace(
107 PDELab::EntitySetPartitioner::Simple{ entity_set, part_patches, halo },
108 TypeTree::treePath()));
109 } else if (part_coloring == "DSatur") {
110 return make_one_step_op(
111 execution_policy,
112 basis.subSpace(PDELab::EntitySetPartitioner::SimpleColored{ entity_set, part_patches, halo },
113 TypeTree::treePath()));
114 } else {
115 throw format_exception(IOError{}, "Not known coloring algorithm '{}' known", part_coloring);
116 }
117 } else {
118 throw format_exception(IOError{}, "Not known partition algorithm '{}' known", part_type);
119 }
120 };
121
122 auto exec_policy = assembly_cfg.get("type", "concurrent");
123#else
124 auto exec_policy = assembly_cfg.get("type", "sequential");
125#endif // DUNE_COPASI_CONCURRENT_ASSEMBLY
126
127 spdlog::info("Creating time-step operator with '{}' execution policy", exec_policy);
128 if (exec_policy != "sequential" and part_patches < std::thread::hardware_concurrency()) {
129 exec_policy = "sequential";
130 spdlog::warn("Too few patches '{}', execution policy of assembler is switched to 'sequential'", part_patches);
131 }
132
133 if (exec_policy == "sequential") {
134 one_step = make_one_step_op(PDELab::Execution::seq, basis);
135#if DUNE_COPASI_CONCURRENT_ASSEMBLY
136 } else if (exec_policy == "concurrent") {
137 one_step = make_concurrent_one_step_op(PDELab::Execution::par);
138#endif // DUNE_COPASI_CONCURRENT_ASSEMBLY
139 } else {
140 throw format_exception(IOError{}, "Not known execution '{}' policy known", exec_policy);
141 }
142
143 return one_step;
144}
145
146} // namespace Dune::Copasi::DiffusionReaction
147
148#endif // DUNE_COPASI_MODEL_MAKE_DIFFUSION_REACTION_STEP_OPERATOR_HH
Definition: factory.hh:28
auto format_exception(Exception &&e, fmt::format_string< Args... > format, Args &&... args)
Definition: exceptions.hh:23