ryujin 2.1.1 revision ef0fcd4010d109b860652ece4a7b8963fb7d46b1
selected_components_extractor.h
Go to the documentation of this file.
1//
2// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
3// Copyright (C) 2024 by the ryujin authors
4//
5
6#pragma once
7
8#include <compile_time_options.h>
9
10#include "state_vector.h"
11
12namespace ryujin
13{
14 template <typename Description, int dim, typename Number>
17
18 using View =
19 typename Description::template HyperbolicSystemView<dim, Number>;
20
22 using InitialPrecomputedVector = typename View::InitialPrecomputedVector;
24
26
27 static void check(const std::vector<std::string> &selected)
28 {
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);
32 };
33
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) ||
39 (entry == "alpha");
40 AssertThrow(
41 found,
42 dealii::ExcMessage(
43 "Invalid component name: \"" + entry +
44 "\" is not a valid conserved/primitive/precomputed/initial "
45 "component name."));
46 }
47 }
48
49 static std::vector<ScalarVector>
50 extract(const HyperbolicSystem &hyperbolic_system,
51 const StateVector &state_vector,
52 const InitialPrecomputedVector &initial_precomputed,
53 const ScalarVector &alpha,
54 const std::vector<std::string> &selected)
55 {
56 /*
57 * Match the selected_components strings against conserved,
58 * primitive, and initial component names and record an index pair
59 * matching return vector position and component index:
60 */
61
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;
67
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});
74 return true;
75 }
76 return false;
77 };
78
79 if (search(View::component_names, conserved_indices))
80 ;
81 else if (search(View::primitive_component_names, primitive_indices))
82 ;
83 else if (search(View::precomputed_names, precomputed_indices))
84 ;
85 else if (search(View::initial_precomputed_names, initial_indices))
86 ;
87 else if (entry == "alpha")
88 alpha_indices.push_back(i++);
89 else
90 AssertThrow(false, dealii::ExcInternalError());
91 }
92
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);
97
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);
101 }
102
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];
112 }
113 }
114
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);
118 }
119
120 for (const auto &[i, k] : initial_indices) {
121 initial_precomputed.extract_component(extracted_components[i], k);
122 }
123
124 for (const auto &i : alpha_indices)
125 extracted_components[i] = alpha;
126
127 return extracted_components;
128 }
129 };
130} // namespace ryujin
dealii::LinearAlgebra::distributed::Vector< Number > ScalarVector
Definition: state_vector.h:33
std::tuple< MultiComponentVector< Number, problem_dim >, MultiComponentVector< Number, prec_dim >, BlockVector< Number > > StateVector
Definition: state_vector.h:53
Euler::HyperbolicSystem HyperbolicSystem
Definition: description.h:34
typename Description::template HyperbolicSystemView< dim, Number > View
typename Description::HyperbolicSystem HyperbolicSystem
Vectors::ScalarVector< Number > ScalarVector
typename View::InitialPrecomputedVector InitialPrecomputedVector
static std::vector< ScalarVector > extract(const HyperbolicSystem &hyperbolic_system, const StateVector &state_vector, const InitialPrecomputedVector &initial_precomputed, const ScalarVector &alpha, const std::vector< std::string > &selected)
static void check(const std::vector< std::string > &selected)