Shamrock 2025.10.0
Astrophysical Code
Loading...
Searching...
No Matches
SlopeLimitedGradient.cpp
Go to the documentation of this file.
1// -------------------------------------------------------//
2//
3// SHAMROCK code for hydrodynamics
4// Copyright (c) 2021-2026 Timothée David--Cléris <tim.shamrock@proton.me>
5// SPDX-License-Identifier: CeCILL Free Software License Agreement v2.1
6// Shamrock is licensed under the CeCILL 2.1 License, see LICENSE for more information
7//
8// -------------------------------------------------------//
9
17
18#include "shammodels/ramses/modules/SlopeLimitedGradient.hpp"
19#include "shammodels/common/amr/NeighGraph.hpp"
20#include "shammodels/ramses/modules/SlopeLimitedGradientUtilities.hpp"
21#include "shamrock/patch/PatchDataField.hpp"
22#include <type_traits>
23
24namespace {
25
26 template<class Tvec, class TgridVec, SlopeMode mode>
27 class KernelSlopeLimScalGrad {
28
29 using Edges = typename shammodels::basegodunov::modules::
30 SlopeLimitedScalarGradient<Tvec, TgridVec>::Edges;
31 using Tscal = shambase::VecComponent<Tvec>;
32 using OrientedAMRGraph = shammodels::basegodunov::modules::OrientedAMRGraph<Tvec, TgridVec>;
33 using AMRGraph = shammodels::basegodunov::modules::AMRGraph;
34
35 public:
36 inline static void kernel(Edges &edges, u32 block_size, u32 var_per_cell) {
37
38 edges.cell_neigh_graph.graph.for_each(
39 [&](u64 id, const OrientedAMRGraph &oriented_cell_graph) {
40 auto &field_span = edges.span_field.get_spans().get(id);
41 auto &field_grad_span = edges.span_grad_field.get_spans().get(id);
42 auto &cell_sizes_span = edges.spans_block_cell_sizes.get_spans().get(id);
43
44 AMRGraph &graph_neigh_xp
45 = shambase::get_check_ref(oriented_cell_graph.graph_links[Direction::xp]);
46 AMRGraph &graph_neigh_xm
47 = shambase::get_check_ref(oriented_cell_graph.graph_links[Direction::xm]);
48 AMRGraph &graph_neigh_yp
49 = shambase::get_check_ref(oriented_cell_graph.graph_links[Direction::yp]);
50 AMRGraph &graph_neigh_ym
51 = shambase::get_check_ref(oriented_cell_graph.graph_links[Direction::ym]);
52 AMRGraph &graph_neigh_zp
53 = shambase::get_check_ref(oriented_cell_graph.graph_links[Direction::zp]);
54 AMRGraph &graph_neigh_zm
55 = shambase::get_check_ref(oriented_cell_graph.graph_links[Direction::zm]);
56
57 sham::EventList depends_list;
58
59 auto cell_sizes = cell_sizes_span.get_read_access(depends_list);
60 auto field = field_span.get_read_access(depends_list);
61 auto field_grad = field_grad_span.get_write_access(depends_list);
62
63 auto graph_iter_xp = graph_neigh_xp.get_read_access(depends_list);
64 auto graph_iter_xm = graph_neigh_xm.get_read_access(depends_list);
65 auto graph_iter_yp = graph_neigh_yp.get_read_access(depends_list);
66 auto graph_iter_ym = graph_neigh_ym.get_read_access(depends_list);
67 auto graph_iter_zp = graph_neigh_zp.get_read_access(depends_list);
68 auto graph_iter_zm = graph_neigh_zm.get_read_access(depends_list);
69
70 sham::DeviceQueue &q = shamsys::instance::get_compute_scheduler().get_queue();
71 auto e = q.submit(depends_list, [&](sycl::handler &cgh) {
72 u32 cell_count = (edges.sizes.indexes.get(id)) * block_size;
73
74 shambase::parallel_for(
75 cgh, cell_count * var_per_cell, "compute_grad_rho", [=](u64 gid) {
76 const u32 tmp_gid = (u32) gid;
77
78 const u32 cell_global_id = tmp_gid / var_per_cell;
79 const u32 var_off_loc = tmp_gid % var_per_cell;
80
81 const u32 block_id = cell_global_id / block_size;
82 const u32 cell_loc_id = cell_global_id % block_size;
83
84 Tscal delta_cell = cell_sizes[block_id];
85
86 auto result = get_3d_grad<Tscal, Tvec, mode>(
87 cell_sizes,
88 block_size,
89 cell_global_id,
90 graph_iter_xp,
91 graph_iter_xm,
92 graph_iter_yp,
93 graph_iter_ym,
94 graph_iter_zp,
95 graph_iter_zm,
96 [=](u32 id) {
97 return field[var_per_cell * id + var_off_loc];
98 });
99
100 field_grad[var_per_cell * cell_global_id + var_off_loc]
101 = {result[0], result[1], result[2]};
102 });
103 });
104
105 cell_sizes_span.complete_event_state(e);
106 field_span.complete_event_state(e);
107 field_grad_span.complete_event_state(e);
108
109 graph_neigh_xp.complete_event_state(e);
110 graph_neigh_xm.complete_event_state(e);
111 graph_neigh_yp.complete_event_state(e);
112 graph_neigh_ym.complete_event_state(e);
113 graph_neigh_zp.complete_event_state(e);
114 graph_neigh_zm.complete_event_state(e);
115 });
116 }
117 };
118
119 template<class Tvec, class TgridVec, SlopeMode mode>
120 class KernelSlopeLimVecGrad {
121
122 using Edges = typename shammodels::basegodunov::modules::
123 SlopeLimitedVectorGradient<Tvec, TgridVec>::Edges;
124 using Tscal = shambase::VecComponent<Tvec>;
125 using OrientedAMRGraph = shammodels::basegodunov::modules::OrientedAMRGraph<Tvec, TgridVec>;
126 using AMRGraph = shammodels::basegodunov::modules::AMRGraph;
127
128 public:
129 inline static void kernel(Edges &edges, u32 block_size, u32 var_per_cell) {
130
131 edges.cell_neigh_graph.graph.for_each(
132 [&](u64 id, const OrientedAMRGraph &oriented_cell_graph) {
133 auto &field_span = edges.span_field.get_spans().get(id);
134 auto &field_dx_span = edges.span_dx_field.get_spans().get(id);
135 auto &field_dy_span = edges.span_dy_field.get_spans().get(id);
136 auto &field_dz_span = edges.span_dz_field.get_spans().get(id);
137 auto &cell_sizes_span = edges.spans_block_cell_sizes.get_spans().get(id);
138
139 AMRGraph &graph_neigh_xp
140 = shambase::get_check_ref(oriented_cell_graph.graph_links[Direction::xp]);
141 AMRGraph &graph_neigh_xm
142 = shambase::get_check_ref(oriented_cell_graph.graph_links[Direction::xm]);
143 AMRGraph &graph_neigh_yp
144 = shambase::get_check_ref(oriented_cell_graph.graph_links[Direction::yp]);
145 AMRGraph &graph_neigh_ym
146 = shambase::get_check_ref(oriented_cell_graph.graph_links[Direction::ym]);
147 AMRGraph &graph_neigh_zp
148 = shambase::get_check_ref(oriented_cell_graph.graph_links[Direction::zp]);
149 AMRGraph &graph_neigh_zm
150 = shambase::get_check_ref(oriented_cell_graph.graph_links[Direction::zm]);
151
152 sham::EventList depends_list;
153
154 auto cell_sizes = cell_sizes_span.get_read_access(depends_list);
155 auto field = field_span.get_read_access(depends_list);
156 auto field_dx = field_dx_span.get_write_access(depends_list);
157 auto field_dy = field_dy_span.get_write_access(depends_list);
158 auto field_dz = field_dz_span.get_write_access(depends_list);
159
160 auto graph_iter_xp = graph_neigh_xp.get_read_access(depends_list);
161 auto graph_iter_xm = graph_neigh_xm.get_read_access(depends_list);
162 auto graph_iter_yp = graph_neigh_yp.get_read_access(depends_list);
163 auto graph_iter_ym = graph_neigh_ym.get_read_access(depends_list);
164 auto graph_iter_zp = graph_neigh_zp.get_read_access(depends_list);
165 auto graph_iter_zm = graph_neigh_zm.get_read_access(depends_list);
166
167 sham::DeviceQueue &q = shamsys::instance::get_compute_scheduler().get_queue();
168 auto e = q.submit(depends_list, [&](sycl::handler &cgh) {
169 u32 cell_count = (edges.sizes.indexes.get(id)) * block_size;
170
171 shambase::parallel_for(
172 cgh, cell_count * var_per_cell, "compute_grad_rho", [=](u64 gid) {
173 const u32 tmp_gid = (u32) gid;
174
175 const u32 cell_global_id = tmp_gid / var_per_cell;
176 const u32 var_off_loc = tmp_gid % var_per_cell;
177
178 const u32 block_id = cell_global_id / block_size;
179 const u32 cell_loc_id = cell_global_id % block_size;
180
181 Tscal delta_cell = cell_sizes[block_id];
182
183 auto result = get_3d_grad<Tvec, Tvec, mode>(
184 cell_sizes,
185 block_size,
186 cell_global_id,
187 graph_iter_xp,
188 graph_iter_xm,
189 graph_iter_yp,
190 graph_iter_ym,
191 graph_iter_zp,
192 graph_iter_zm,
193 [=](u32 id) {
194 return field[var_per_cell * id + var_off_loc];
195 });
196
197 field_dx[var_per_cell * cell_global_id + var_off_loc] = result[0];
198 field_dy[var_per_cell * cell_global_id + var_off_loc] = result[1];
199 field_dz[var_per_cell * cell_global_id + var_off_loc] = result[2];
200 });
201 });
202
203 cell_sizes_span.complete_event_state(e);
204 field_span.complete_event_state(e);
205 field_dx_span.complete_event_state(e);
206 field_dy_span.complete_event_state(e);
207 field_dz_span.complete_event_state(e);
208
209 graph_neigh_xp.complete_event_state(e);
210 graph_neigh_xm.complete_event_state(e);
211 graph_neigh_yp.complete_event_state(e);
212 graph_neigh_ym.complete_event_state(e);
213 graph_neigh_zp.complete_event_state(e);
214 graph_neigh_zm.complete_event_state(e);
215 });
216 }
217 };
218} // namespace
219
220namespace shammodels::basegodunov::modules {
221
222 template<class Tvec, class TgridVec>
223 void SlopeLimitedScalarGradient<Tvec, TgridVec>::_impl_evaluate_internal() {
224 StackEntry stack_loc{};
225 auto edges = get_edges();
226
227 edges.spans_block_cell_sizes.check_sizes(edges.sizes.indexes);
228 edges.span_field.check_sizes(edges.sizes.indexes);
229
230 edges.span_grad_field.ensure_sizes(edges.sizes.indexes);
231
232 if (mode == SlopeMode::None) {
233 using Kern = KernelSlopeLimScalGrad<Tvec, TgridVec, None>;
234 Kern::kernel(edges, block_size, var_per_cell);
235 } else if (mode == SlopeMode::VanLeer_f) {
236 using Kern = KernelSlopeLimScalGrad<Tvec, TgridVec, VanLeer_f>;
237 Kern::kernel(edges, block_size, var_per_cell);
238 } else if (mode == SlopeMode::VanLeer_std) {
239 using Kern = KernelSlopeLimScalGrad<Tvec, TgridVec, VanLeer_std>;
240 Kern::kernel(edges, block_size, var_per_cell);
241 } else if (mode == SlopeMode::VanLeer_sym) {
242 using Kern = KernelSlopeLimScalGrad<Tvec, TgridVec, VanLeer_sym>;
243 Kern::kernel(edges, block_size, var_per_cell);
244 } else if (mode == SlopeMode::Minmod) {
245 using Kern = KernelSlopeLimScalGrad<Tvec, TgridVec, Minmod>;
246 Kern::kernel(edges, block_size, var_per_cell);
247 } else {
248 shambase::throw_unimplemented();
249 }
250 }
251
252 template<class Tvec, class TgridVec>
253 std::string SlopeLimitedScalarGradient<Tvec, TgridVec>::_impl_get_tex() const {
254
255 std::string sizes = get_ro_edge_base(0).get_tex_symbol();
256 std::string cell_neigh_graph = get_ro_edge_base(1).get_tex_symbol();
257 std::string spans_block_cell_sizes = get_ro_edge_base(2).get_tex_symbol();
258 std::string span_field = get_ro_edge_base(3).get_tex_symbol();
259 std::string span_grad_field = get_rw_edge_base(0).get_tex_symbol();
260
261 std::string tex = R"tex(
262 Slope limited gradient (Scalar)
263 )tex";
264
265 shambase::replace_all(tex, "{sizes}", sizes);
266 shambase::replace_all(tex, "{cell_neigh_graph}", cell_neigh_graph);
267 shambase::replace_all(tex, "{spans_block_cell_sizes}", spans_block_cell_sizes);
268 shambase::replace_all(tex, "{span_field}", span_field);
269 shambase::replace_all(tex, "{span_grad_field}", span_grad_field);
270
271 return tex;
272 }
273 template<class Tvec, class TgridVec>
274 void SlopeLimitedVectorGradient<Tvec, TgridVec>::_impl_evaluate_internal() {
275 StackEntry stack_loc{};
276 auto edges = get_edges();
277
278 edges.spans_block_cell_sizes.check_sizes(edges.sizes.indexes);
279 edges.span_field.check_sizes(edges.sizes.indexes);
280
281 edges.span_dx_field.ensure_sizes(edges.sizes.indexes);
282 edges.span_dy_field.ensure_sizes(edges.sizes.indexes);
283 edges.span_dz_field.ensure_sizes(edges.sizes.indexes);
284
285 if (mode == SlopeMode::None) {
286 using Kern = KernelSlopeLimVecGrad<Tvec, TgridVec, None>;
287 Kern::kernel(edges, block_size, var_per_cell);
288 } else if (mode == SlopeMode::VanLeer_f) {
289 using Kern = KernelSlopeLimVecGrad<Tvec, TgridVec, VanLeer_f>;
290 Kern::kernel(edges, block_size, var_per_cell);
291 } else if (mode == SlopeMode::VanLeer_std) {
292 using Kern = KernelSlopeLimVecGrad<Tvec, TgridVec, VanLeer_std>;
293 Kern::kernel(edges, block_size, var_per_cell);
294 } else if (mode == SlopeMode::VanLeer_sym) {
295 using Kern = KernelSlopeLimVecGrad<Tvec, TgridVec, VanLeer_sym>;
296 Kern::kernel(edges, block_size, var_per_cell);
297 } else if (mode == SlopeMode::Minmod) {
298 using Kern = KernelSlopeLimVecGrad<Tvec, TgridVec, Minmod>;
299 Kern::kernel(edges, block_size, var_per_cell);
300 } else {
301 shambase::throw_unimplemented();
302 }
303 }
304
305 template<class Tvec, class TgridVec>
306 std::string SlopeLimitedVectorGradient<Tvec, TgridVec>::_impl_get_tex() const {
307
308 std::string sizes = get_ro_edge_base(0).get_tex_symbol();
309 std::string cell_neigh_graph = get_ro_edge_base(1).get_tex_symbol();
310 std::string spans_block_cell_sizes = get_ro_edge_base(2).get_tex_symbol();
311 std::string span_field = get_ro_edge_base(3).get_tex_symbol();
312 std::string span_dx_field = get_rw_edge_base(0).get_tex_symbol();
313 std::string span_dy_field = get_rw_edge_base(0).get_tex_symbol();
314 std::string span_dz_field = get_rw_edge_base(0).get_tex_symbol();
315
316 std::string tex = R"tex(
317 Slope limited gradient (Vector)
318 )tex";
319
320 shambase::replace_all(tex, "{sizes}", sizes);
321 shambase::replace_all(tex, "{cell_neigh_graph}", cell_neigh_graph);
322 shambase::replace_all(tex, "{spans_block_cell_sizes}", spans_block_cell_sizes);
323 shambase::replace_all(tex, "{span_field}", span_field);
324 shambase::replace_all(tex, "{span_dx_field}", span_dx_field);
325 shambase::replace_all(tex, "{span_dy_field}", span_dy_field);
326 shambase::replace_all(tex, "{span_dz_field}", span_dz_field);
327
328 return tex;
329 }
330
331} // namespace shammodels::basegodunov::modules
332
333template class shammodels::basegodunov::modules::SlopeLimitedScalarGradient<f64_3, i64_3>;
334template class shammodels::basegodunov::modules::SlopeLimitedVectorGradient<f64_3, i64_3>;
u32
std::uint32_t u32
32 bit unsigned integer
Definition aliases_int.hpp:27
u64
std::uint64_t u64
64 bit unsigned integer
Definition aliases_int.hpp:26
sham::DeviceQueue::submit
sycl::event submit(Fct &&fct)
Submits a kernel to the SYCL queue.
Definition DeviceQueue.hpp:101
sham::DeviceScheduler::get_queue
DeviceQueue & get_queue(u32 id=0)
Get a reference to a DeviceQueue.
Definition DeviceScheduler.cpp:57
shammodels::basegodunov::modules::SlopeLimitedScalarGradient::_impl_get_tex
virtual std::string _impl_get_tex() const
get the tex of the node
Definition SlopeLimitedGradient.cpp:253
shammodels::basegodunov::modules::SlopeLimitedVectorGradient::_impl_get_tex
virtual std::string _impl_get_tex() const
get the tex of the node
Definition SlopeLimitedGradient.cpp:304
shamrock::solvergraph::INode::get_rw_edge_base
IEdge & get_rw_edge_base(int slot)
Get a reference to a read write edge and cast it to the type IEdge.
Definition INode.hpp:137
shamrock::solvergraph::INode::get_ro_edge_base
const IEdge & get_ro_edge_base(int slot)
Get a reference to a read only edge.
Definition INode.hpp:128
shambase::replace_all
void replace_all(std::string &inout, std::string_view what, std::string_view with)
replace all occurence of a search string with another
Definition string.hpp:110
shambase::get_check_ref
T & get_check_ref(const std::unique_ptr< T > &ptr, SourceLocation loc=SourceLocation())
Takes a std::unique_ptr and returns a reference to the object it holds. It throws a std::runtime_erro...
Definition memory.hpp:110
shambase::throw_unimplemented
void throw_unimplemented(SourceLocation loc=SourceLocation{})
Throw a std::runtime_error saying that the function is unimplemented.
Definition exception.hpp:153
shammodels::basegodunov::modules
namespace for the basegodunov model modules
Definition AMRCellInfos.hpp:20
StackEntry
shambase::details::BasicStackEntry StackEntry
Alias for shambase::details::BasicStackEntry.
Definition stacktrace.hpp:242
Generated on for Shamrock by doxygen 1.15.0