8#include <compile_time_options.h>
14 template <
typename Description,
int dim,
typename Number>
27 static void check(
const std::vector<std::string> &selected)
29 const auto search = [&](
const auto entry,
const auto &names) {
30 const auto pos = std::find(std::begin(names), std::end(names), entry);
31 return pos != std::end(names);
34 for (
const auto &entry : selected) {
35 const auto found = search(entry, View::component_names) ||
36 search(entry, View::primitive_component_names) ||
37 search(entry, View::precomputed_names) ||
38 search(entry, View::initial_precomputed_names) ||
43 "Invalid component name: \"" + entry +
44 "\" is not a valid conserved/primitive/precomputed/initial "
49 static std::vector<ScalarVector>
54 const std::vector<std::string> &selected)
62 std::vector<std::tuple<std::size_t, std::size_t>> conserved_indices;
63 std::vector<std::tuple<std::size_t, std::size_t>> primitive_indices;
64 std::vector<std::tuple<std::size_t, std::size_t>> precomputed_indices;
65 std::vector<std::tuple<std::size_t, std::size_t>> initial_indices;
66 std::vector<std::size_t> alpha_indices;
68 for (std::size_t i = 0;
const auto &entry : selected) {
69 const auto search = [&](
const auto &names,
auto &indices) {
70 const auto pos = std::find(std::begin(names), std::end(names), entry);
71 if (pos != std::end(names)) {
72 const auto index = std::distance(std::begin(names), pos);
73 indices.push_back({i++, index});
79 if (search(View::component_names, conserved_indices))
81 else if (search(View::primitive_component_names, primitive_indices))
83 else if (search(View::precomputed_names, precomputed_indices))
85 else if (search(View::initial_precomputed_names, initial_indices))
87 else if (entry ==
"alpha")
88 alpha_indices.push_back(i++);
90 AssertThrow(
false, dealii::ExcInternalError());
93 std::vector<ScalarVector> extracted_components(selected.size());
94 const auto &scalar_partitioner = alpha.get_partitioner();
95 for (
auto &it : extracted_components)
96 it.reinit(scalar_partitioner);
98 for (
const auto &[i, k] : conserved_indices) {
99 const auto &U = std::get<0>(state_vector);
100 U.extract_component(extracted_components[i], k);
103 if (!primitive_indices.empty()) {
104 const auto &U = std::get<0>(state_vector);
105 const unsigned int n_owned = scalar_partitioner->locally_owned_size();
106 const auto view = hyperbolic_system.template view<dim, Number>();
107 for (
unsigned int i = 0; i < n_owned; ++i) {
108 const auto U_i = U.get_tensor(i);
109 const auto PU_i = view.to_primitive_state(U_i);
110 for (
const auto &[j, k] : primitive_indices)
111 extracted_components[j].local_element(i) = PU_i[k];
115 for (
const auto &[i, k] : precomputed_indices) {
116 const auto &prec = std::get<1>(state_vector);
117 prec.extract_component(extracted_components[i], k);
120 for (
const auto &[i, k] : initial_indices) {
121 initial_precomputed.extract_component(extracted_components[i], k);
124 for (
const auto &i : alpha_indices)
125 extracted_components[i] = alpha;
127 return extracted_components;
dealii::LinearAlgebra::distributed::Vector< Number > ScalarVector
std::tuple< MultiComponentVector< Number, problem_dim >, MultiComponentVector< Number, prec_dim >, BlockVector< Number > > StateVector
Euler::HyperbolicSystem HyperbolicSystem