Shamrock 2025.10.0
Astrophysical Code
Loading...
Searching...
No Matches
UpdateDerivs.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
29#include "shammodels/gsph/modules/UpdateDerivs.hpp"
30#include "shambackends/math.hpp"
31#include "shammath/sphkernels.hpp"
32#include "shammodels/gsph/config/FieldNames.hpp"
33#include "shammodels/gsph/math/forces.hpp"
34#include "shammodels/gsph/math/riemann/iterative.hpp"
35#include "shammodels/sph/math/density.hpp"
36#include "shamsys/legacy/log.hpp"
37#include "shamtree/TreeTraversal.hpp"
38
39template<class Tvec, template<class> class SPHKernel>
40void shammodels::gsph::modules::UpdateDerivs<Tvec, SPHKernel>::update_derivs() {
41 StackEntry stack_loc{};
42
43 Cfg_Riemann cfg_riemann = solver_config.riemann_config;
44
45 if (Iterative *v = std::get_if<Iterative>(&cfg_riemann.config)) {
46 update_derivs_iterative(*v);
47 } else if (HLLC *v = std::get_if<HLLC>(&cfg_riemann.config)) {
48 update_derivs_hllc(*v);
49 } else {
50 shambase::throw_unimplemented("Riemann solver type not supported by UpdateDerivs");
51 }
52}
53
54template<class Tvec, template<class> class SPHKernel>
55void shammodels::gsph::modules::UpdateDerivs<Tvec, SPHKernel>::update_derivs_iterative(
56 Iterative cfg) {
57
58 StackEntry stack_loc{};
59
60 using namespace shamrock;
61 using namespace shamrock::patch;
62
63 PatchDataLayerLayout &pdl = scheduler().pdl_old();
64
65 // Get field indices from the patch data layout
66 const u32 ixyz = pdl.get_field_idx<Tvec>(gsph::names::common::xyz);
67 const u32 ivxyz = pdl.get_field_idx<Tvec>(gsph::names::newtonian::vxyz);
68 const u32 iaxyz = pdl.get_field_idx<Tvec>(gsph::names::newtonian::axyz);
69 const u32 ihpart = pdl.get_field_idx<Tscal>(gsph::names::common::hpart);
70
71 // Optional internal energy fields (for adiabatic EOS)
72 const bool has_uint = solver_config.has_field_uint();
73 const u32 iuint = has_uint ? pdl.get_field_idx<Tscal>(gsph::names::newtonian::uint) : 0;
74 const u32 iduint = has_uint ? pdl.get_field_idx<Tscal>(gsph::names::newtonian::duint) : 0;
75
76 // Ghost layout for neighbor data
77 shamrock::patch::PatchDataLayerLayout &ghost_layout
78 = shambase::get_check_ref(storage.ghost_layout.get());
79 u32 ihpart_interf = ghost_layout.get_field_idx<Tscal>(gsph::names::common::hpart);
80 u32 ivxyz_interf = ghost_layout.get_field_idx<Tvec>(gsph::names::newtonian::vxyz);
81 u32 iomega_interf = ghost_layout.get_field_idx<Tscal>(gsph::names::newtonian::omega);
82 u32 idensity_interf = ghost_layout.get_field_idx<Tscal>(gsph::names::newtonian::density);
83 u32 iuint_interf
84 = has_uint ? ghost_layout.get_field_idx<Tscal>(gsph::names::newtonian::uint) : 0;
85
86 // Get merged data and caches from storage
87 auto &merged_xyzh = storage.merged_xyzh.get();
88 shamrock::solvergraph::Field<Tscal> &omega_field = shambase::get_check_ref(storage.omega);
89 shambase::DistributedData<PatchDataLayer> &mpdats = storage.merged_patchdata_ghost.get();
90
91 // Get pressure and soundspeed from storage (includes ghosts)
92 shamrock::solvergraph::Field<Tscal> &pressure_field = shambase::get_check_ref(storage.pressure);
93 shamrock::solvergraph::Field<Tscal> &soundspeed_field
94 = shambase::get_check_ref(storage.soundspeed);
95
96 // Iterate over all non-empty patches
97 scheduler().for_each_patchdata_nonempty([&](Patch cur_p, PatchDataLayer &pdat) {
98 PatchDataLayer &mpdat = mpdats.get(cur_p.id_patch);
99
100 // Get buffers for local and ghost data
101 sham::DeviceBuffer<Tvec> &buf_xyz
102 = merged_xyzh.get(cur_p.id_patch).template get_field_buf_ref<Tvec>(0);
103 sham::DeviceBuffer<Tvec> &buf_axyz = pdat.get_field_buf_ref<Tvec>(iaxyz);
104 sham::DeviceBuffer<Tvec> &buf_vxyz = mpdat.get_field_buf_ref<Tvec>(ivxyz_interf);
105 sham::DeviceBuffer<Tscal> &buf_hpart = mpdat.get_field_buf_ref<Tscal>(ihpart_interf);
106 sham::DeviceBuffer<Tscal> &buf_omega = mpdat.get_field_buf_ref<Tscal>(iomega_interf);
107 sham::DeviceBuffer<Tscal> &buf_pressure
108 = pressure_field.get_field(cur_p.id_patch).get_buf();
109 sham::DeviceBuffer<Tscal> &buf_cs = soundspeed_field.get_field(cur_p.id_patch).get_buf();
110
111 // Get neighbor cache for this patch
112 tree::ObjectCache &pcache
113 = shambase::get_check_ref(storage.neigh_cache).get_cache(cur_p.id_patch);
114
115 // Set up SYCL queue and event tracking
116 sham::DeviceQueue &q = shamsys::instance::get_compute_scheduler().get_queue();
117 sham::EventList depends_list;
118
119 // Get density from merged ghost data (SPH summation density)
120 sham::DeviceBuffer<Tscal> &buf_density = mpdat.get_field_buf_ref<Tscal>(idensity_interf);
121
122 // Get buffer accessors
123 auto xyz = buf_xyz.get_read_access(depends_list);
124 auto axyz = buf_axyz.get_write_access(depends_list);
125 auto vxyz = buf_vxyz.get_read_access(depends_list);
126 auto hpart = buf_hpart.get_read_access(depends_list);
127 auto omega_acc = buf_omega.get_read_access(depends_list);
128 auto density_acc = buf_density.get_read_access(depends_list);
129 auto pressure_acc = buf_pressure.get_read_access(depends_list);
130 auto cs_acc = buf_cs.get_read_access(depends_list);
131 auto ploop_ptrs = pcache.get_read_access(depends_list);
132
133 // Optional: internal energy
134 sham::DeviceBuffer<Tscal> *buf_duint_ptr = nullptr;
135 Tscal *duint_acc = nullptr;
136 if (has_uint) {
137 buf_duint_ptr = &pdat.get_field_buf_ref<Tscal>(iduint);
138 duint_acc = buf_duint_ptr->get_write_access(depends_list);
139 }
140
141 auto e = q.submit(depends_list, [&](sycl::handler &cgh) {
142 const Tscal pmass = solver_config.gpart_mass;
143 const Tscal gamma = solver_config.get_eos_gamma();
144 const Tscal tol = cfg.tol;
145 const u32 max_iter = cfg.max_iter;
146 const bool do_energy = has_uint;
147
148 // Use shamrock's ObjectCacheIterator for neighbor traversal
149 tree::ObjectCacheIterator particle_looper(ploop_ptrs);
150
151 constexpr Tscal Rker2 = Kernel::Rkern * Kernel::Rkern;
152
153 shambase::parallel_for(cgh, pdat.get_obj_cnt(), "GSPH derivs iterative", [=](u64 gid) {
154 u32 id_a = (u32) gid;
155
156 using namespace shamrock::sph;
157
158 // Initialize accumulators
159 Tvec sum_axyz = {0, 0, 0};
160 Tscal sum_du_a = 0;
161
162 // Particle a state
163 const Tscal h_a = hpart[id_a];
164 const Tvec xyz_a = xyz[id_a];
165 const Tvec vxyz_a = vxyz[id_a];
166 const Tscal omega_a = omega_acc[id_a];
167
168 // Use SPH-summation density
169 const Tscal rho_a = sycl::max(density_acc[id_a], Tscal(1e-30));
170
171 // Use pressure and soundspeed from storage
172 const Tscal P_a = sycl::max(pressure_acc[id_a], Tscal(1e-30));
173 const Tscal cs_a = sycl::max(cs_acc[id_a], Tscal(1e-10));
174
175 // Loop over neighbors
176 particle_looper.for_each_object(id_a, [&](u32 id_b) {
177 if (id_a == id_b)
178 return; // Skip self
179
180 // Distance and kernel support check
181 const Tvec dr = xyz_a - xyz[id_b];
182 const Tscal rab2 = sycl::dot(dr, dr);
183 const Tscal h_b = hpart[id_b];
184
185 // Skip if outside kernel support
186 if (rab2 > h_a * h_a * Rker2 && rab2 > h_b * h_b * Rker2) {
187 return;
188 }
189
190 const Tscal rab = sycl::sqrt(rab2);
191 const Tvec vxyz_b = vxyz[id_b];
192 const Tscal omega_b = omega_acc[id_b];
193
194 // Use SPH-summation density
195 const Tscal rho_b = sycl::max(density_acc[id_b], Tscal(1e-30));
196
197 // Use pressure and soundspeed from storage
198 const Tscal P_b = sycl::max(pressure_acc[id_b], Tscal(1e-30));
199 const Tscal cs_b = sycl::max(cs_acc[id_b], Tscal(1e-10));
200
201 // Unit vector from a to b
202 const Tscal rab_inv = sham::inv_sat_positive(rab);
203 const Tvec r_ab_unit = dr * rab_inv;
204
205 // Project velocities onto pair axis for 1D Riemann problem
206 const Tscal u_a_proj = sycl::dot(vxyz_a, r_ab_unit);
207 const Tscal u_b_proj = sycl::dot(vxyz_b, r_ab_unit);
208
209 // Solve 1D Riemann problem using iterative solver
210 // Convention: Left state = neighbor b, Right state = current a
211 auto riemann_result = riemann::iterative_solver<Tscal>(
212 u_b_proj,
213 rho_b,
214 P_b, // Left = neighbor
215 u_a_proj,
216 rho_a,
217 P_a, // Right = current
218 gamma,
219 tol,
220 max_iter);
221 const Tscal p_star = riemann_result.p_star;
222 const Tscal v_star = riemann_result.v_star;
223
224 // Kernel gradients
225 const Tscal Fab_a = Kernel::dW_3d(rab, h_a);
226 const Tscal Fab_b = Kernel::dW_3d(rab, h_b);
227
228 // GSPH force contribution
229 shammodels::gsph::add_gsph_force_contribution<Tvec, Tscal>(
230 pmass,
231 p_star,
232 v_star,
233 rho_a,
234 rho_b,
235 omega_a,
236 omega_b,
237 Fab_a,
238 Fab_b,
239 r_ab_unit,
240 vxyz_a,
241 sum_axyz,
242 sum_du_a);
243 });
244
245 // Write accumulated derivatives
246 axyz[id_a] = sum_axyz;
247 if (duint_acc != nullptr) {
248 duint_acc[id_a] = sum_du_a;
249 }
250 });
251 });
252
253 // Complete event states for all buffers
254 buf_xyz.complete_event_state(e);
255 buf_axyz.complete_event_state(e);
256 buf_vxyz.complete_event_state(e);
257 buf_hpart.complete_event_state(e);
258 buf_omega.complete_event_state(e);
259 buf_density.complete_event_state(e);
260 buf_pressure.complete_event_state(e);
261 buf_cs.complete_event_state(e);
262
263 if (has_uint && buf_duint_ptr) {
264 buf_duint_ptr->complete_event_state(e);
265 }
266
267 sham::EventList resulting_events;
268 resulting_events.add_event(e);
269 pcache.complete_event_state(resulting_events);
270 });
271}
272
273template<class Tvec, template<class> class SPHKernel>
274void shammodels::gsph::modules::UpdateDerivs<Tvec, SPHKernel>::update_derivs_hllc(HLLC cfg) {
275
276 StackEntry stack_loc{};
277
278 using namespace shamrock;
279 using namespace shamrock::patch;
280
281 PatchDataLayerLayout &pdl = scheduler().pdl_old();
282
283 // Get field indices
284 const u32 ixyz = pdl.get_field_idx<Tvec>(gsph::names::common::xyz);
285 const u32 ivxyz = pdl.get_field_idx<Tvec>(gsph::names::newtonian::vxyz);
286 const u32 iaxyz = pdl.get_field_idx<Tvec>(gsph::names::newtonian::axyz);
287 const u32 ihpart = pdl.get_field_idx<Tscal>(gsph::names::common::hpart);
288
289 const bool has_uint = solver_config.has_field_uint();
290 const u32 iuint = has_uint ? pdl.get_field_idx<Tscal>(gsph::names::newtonian::uint) : 0;
291 const u32 iduint = has_uint ? pdl.get_field_idx<Tscal>(gsph::names::newtonian::duint) : 0;
292
293 // Ghost layout
294 shamrock::patch::PatchDataLayerLayout &ghost_layout
295 = shambase::get_check_ref(storage.ghost_layout.get());
296 u32 ihpart_interf = ghost_layout.get_field_idx<Tscal>(gsph::names::common::hpart);
297 u32 ivxyz_interf = ghost_layout.get_field_idx<Tvec>(gsph::names::newtonian::vxyz);
298 u32 iomega_interf = ghost_layout.get_field_idx<Tscal>(gsph::names::newtonian::omega);
299 u32 idensity_interf = ghost_layout.get_field_idx<Tscal>(gsph::names::newtonian::density);
300 u32 iuint_interf
301 = has_uint ? ghost_layout.get_field_idx<Tscal>(gsph::names::newtonian::uint) : 0;
302
303 auto &merged_xyzh = storage.merged_xyzh.get();
304 shamrock::solvergraph::Field<Tscal> &omega_field = shambase::get_check_ref(storage.omega);
305 shambase::DistributedData<PatchDataLayer> &mpdats = storage.merged_patchdata_ghost.get();
306
307 // Get pressure and soundspeed from storage (includes ghosts)
308 shamrock::solvergraph::Field<Tscal> &pressure_field = shambase::get_check_ref(storage.pressure);
309 shamrock::solvergraph::Field<Tscal> &soundspeed_field
310 = shambase::get_check_ref(storage.soundspeed);
311
312 scheduler().for_each_patchdata_nonempty([&](Patch cur_p, PatchDataLayer &pdat) {
313 PatchDataLayer &mpdat = mpdats.get(cur_p.id_patch);
314
315 sham::DeviceBuffer<Tvec> &buf_xyz
316 = merged_xyzh.get(cur_p.id_patch).template get_field_buf_ref<Tvec>(0);
317 sham::DeviceBuffer<Tvec> &buf_axyz = pdat.get_field_buf_ref<Tvec>(iaxyz);
318 sham::DeviceBuffer<Tvec> &buf_vxyz = mpdat.get_field_buf_ref<Tvec>(ivxyz_interf);
319 sham::DeviceBuffer<Tscal> &buf_hpart = mpdat.get_field_buf_ref<Tscal>(ihpart_interf);
320 sham::DeviceBuffer<Tscal> &buf_omega = mpdat.get_field_buf_ref<Tscal>(iomega_interf);
321 sham::DeviceBuffer<Tscal> &buf_pressure
322 = pressure_field.get_field(cur_p.id_patch).get_buf();
323 sham::DeviceBuffer<Tscal> &buf_cs = soundspeed_field.get_field(cur_p.id_patch).get_buf();
324
325 tree::ObjectCache &pcache
326 = shambase::get_check_ref(storage.neigh_cache).get_cache(cur_p.id_patch);
327
328 sham::DeviceQueue &q = shamsys::instance::get_compute_scheduler().get_queue();
329 sham::EventList depends_list;
330
331 // Get density from merged ghost data
332 sham::DeviceBuffer<Tscal> &buf_density = mpdat.get_field_buf_ref<Tscal>(idensity_interf);
333
334 auto xyz = buf_xyz.get_read_access(depends_list);
335 auto axyz = buf_axyz.get_write_access(depends_list);
336 auto vxyz = buf_vxyz.get_read_access(depends_list);
337 auto hpart = buf_hpart.get_read_access(depends_list);
338 auto omega_acc = buf_omega.get_read_access(depends_list);
339 auto density_acc = buf_density.get_read_access(depends_list);
340 auto pressure_acc = buf_pressure.get_read_access(depends_list);
341 auto cs_acc = buf_cs.get_read_access(depends_list);
342 auto ploop_ptrs = pcache.get_read_access(depends_list);
343
344 sham::DeviceBuffer<Tscal> *buf_duint_ptr = nullptr;
345 Tscal *duint_acc = nullptr;
346 if (has_uint) {
347 buf_duint_ptr = &pdat.get_field_buf_ref<Tscal>(iduint);
348 duint_acc = buf_duint_ptr->get_write_access(depends_list);
349 }
350
351 auto e = q.submit(depends_list, [&](sycl::handler &cgh) {
352 const Tscal pmass = solver_config.gpart_mass;
353 const Tscal gamma = solver_config.get_eos_gamma();
354 const bool do_energy = has_uint;
355
356 tree::ObjectCacheIterator particle_looper(ploop_ptrs);
357
358 constexpr Tscal Rker2 = Kernel::Rkern * Kernel::Rkern;
359
360 shambase::parallel_for(cgh, pdat.get_obj_cnt(), "GSPH derivs HLLC", [=](u64 gid) {
361 u32 id_a = (u32) gid;
362
363 using namespace shamrock::sph;
364
365 Tvec sum_axyz = {0, 0, 0};
366 Tscal sum_du_a = 0;
367
368 const Tscal h_a = hpart[id_a];
369 const Tvec xyz_a = xyz[id_a];
370 const Tvec vxyz_a = vxyz[id_a];
371 const Tscal omega_a = omega_acc[id_a];
372
373 // Use SPH-summation density
374 const Tscal rho_a = sycl::max(density_acc[id_a], Tscal(1e-30));
375
376 // Use pressure and soundspeed from storage
377 const Tscal P_a = sycl::max(pressure_acc[id_a], Tscal(1e-30));
378 const Tscal cs_a = sycl::max(cs_acc[id_a], Tscal(1e-10));
379
380 particle_looper.for_each_object(id_a, [&](u32 id_b) {
381 if (id_a == id_b)
382 return;
383
384 const Tvec dr = xyz_a - xyz[id_b];
385 const Tscal rab2 = sycl::dot(dr, dr);
386 const Tscal h_b = hpart[id_b];
387
388 if (rab2 > h_a * h_a * Rker2 && rab2 > h_b * h_b * Rker2) {
389 return;
390 }
391
392 const Tscal rab = sycl::sqrt(rab2);
393 const Tvec vxyz_b = vxyz[id_b];
394 const Tscal omega_b = omega_acc[id_b];
395
396 // Use SPH-summation density
397 const Tscal rho_b = sycl::max(density_acc[id_b], Tscal(1e-30));
398
399 // Use pressure and soundspeed from storage
400 const Tscal P_b = sycl::max(pressure_acc[id_b], Tscal(1e-30));
401 const Tscal cs_b = sycl::max(cs_acc[id_b], Tscal(1e-10));
402
403 const Tscal rab_inv = sham::inv_sat_positive(rab);
404 const Tvec r_ab_unit = dr * rab_inv;
405
406 // Project velocities onto pair axis for 1D Riemann problem
407 const Tscal u_a_proj = sycl::dot(vxyz_a, r_ab_unit);
408 const Tscal u_b_proj = sycl::dot(vxyz_b, r_ab_unit);
409
410 // Use HLLC approximate Riemann solver
411 // Convention: Left state = neighbor b, Right state = current a
412 auto riemann_result = riemann::hllc_solver<Tscal>(
413 u_b_proj,
414 rho_b,
415 P_b, // Left = neighbor
416 u_a_proj,
417 rho_a,
418 P_a, // Right = current
419 gamma);
420 const Tscal p_star = riemann_result.p_star;
421 const Tscal v_star = riemann_result.v_star;
422
423 const Tscal Fab_a = Kernel::dW_3d(rab, h_a);
424 const Tscal Fab_b = Kernel::dW_3d(rab, h_b);
425
426 // GSPH force contribution
427 shammodels::gsph::add_gsph_force_contribution<Tvec, Tscal>(
428 pmass,
429 p_star,
430 v_star,
431 rho_a,
432 rho_b,
433 omega_a,
434 omega_b,
435 Fab_a,
436 Fab_b,
437 r_ab_unit,
438 vxyz_a,
439 sum_axyz,
440 sum_du_a);
441 });
442
443 axyz[id_a] = sum_axyz;
444 if (duint_acc != nullptr) {
445 duint_acc[id_a] = sum_du_a;
446 }
447 });
448 });
449
450 buf_xyz.complete_event_state(e);
451 buf_axyz.complete_event_state(e);
452 buf_vxyz.complete_event_state(e);
453 buf_hpart.complete_event_state(e);
454 buf_omega.complete_event_state(e);
455 buf_density.complete_event_state(e);
456 buf_pressure.complete_event_state(e);
457 buf_cs.complete_event_state(e);
458
459 if (has_uint && buf_duint_ptr) {
460 buf_duint_ptr->complete_event_state(e);
461 }
462
463 sham::EventList resulting_events;
464 resulting_events.add_event(e);
465 pcache.complete_event_state(resulting_events);
466 });
467}
468
469// Explicit template instantiations
470// M-spline kernels (Monaghan)
471using namespace shammath;
472template class shammodels::gsph::modules::UpdateDerivs<f64_3, M4>;
473template class shammodels::gsph::modules::UpdateDerivs<f64_3, M6>;
474template class shammodels::gsph::modules::UpdateDerivs<f64_3, M8>;
475
476// Wendland kernels (C2, C4, C6)
477template class shammodels::gsph::modules::UpdateDerivs<f64_3, C2>;
478template class shammodels::gsph::modules::UpdateDerivs<f64_3, C4>;
479template class shammodels::gsph::modules::UpdateDerivs<f64_3, C6>;
FieldNames.hpp
Constants for field names in GSPH solver, organized by physics mode.
shammodels::gsph::names::newtonian::axyz
constexpr const char * axyz
3-acceleration field
Definition FieldNames.hpp:46
shammodels::gsph::names::newtonian::vxyz
constexpr const char * vxyz
3-velocity field
Definition FieldNames.hpp:43
shammodels::gsph::names::common::xyz
constexpr const char * xyz
Position field (3D coordinates)
Definition FieldNames.hpp:33
shammodels::gsph::names::common::hpart
constexpr const char * hpart
Smoothing length field.
Definition FieldNames.hpp:34
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::DeviceBuffer
A buffer allocated in USM (Unified Shared Memory)
Definition DeviceBuffer.hpp:106
sham::DeviceBuffer::complete_event_state
void complete_event_state(sycl::event e) const
Complete the event state of the buffer.
Definition DeviceBuffer.hpp:368
sham::DeviceBuffer::get_write_access
T * get_write_access(sham::EventList &depends_list, SourceLocation src_loc=SourceLocation{})
Get a read-write pointer to the buffer's data.
Definition DeviceBuffer.hpp:349
sham::DeviceBuffer::get_read_access
const T * get_read_access(sham::EventList &depends_list, SourceLocation src_loc=SourceLocation{}) const
Get a read-only pointer to the buffer's data.
Definition DeviceBuffer.hpp:331
sham::DeviceQueue
A SYCL queue associated with a device and a context.
Definition DeviceQueue.hpp:32
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
sham::EventList
Class to manage a list of SYCL events.
Definition EventList.hpp:31
sham::EventList::add_event
void add_event(sycl::event e)
Add an event to the list of events.
Definition EventList.hpp:87
shambase::DistributedData
Represents a collection of objects distributed across patches identified by a u64 id.
Definition DistributedData.hpp:43
shambase::DistributedData::get
T & get(u64 id)
Returns a reference to an object in the collection.
Definition DistributedData.hpp:145
shammodels::gsph::modules::UpdateDerivs
GSPH derivative update module.
Definition UpdateDerivs.hpp:52
shammodels::gsph::modules::UpdateDerivs::update_derivs
void update_derivs()
Update all derivatives using GSPH Riemann solver approach.
Definition UpdateDerivs.cpp:40
shamrock::patch::PatchDataLayerLayout::get_field_idx
u32 get_field_idx(const std::string &field_name) const
Get the field id if matching name & type.
Definition PatchDataLayerLayout.hpp:349
shamrock::patch::PatchDataLayer
PatchDataLayer container class, the layout is described in patchdata_layout.
Definition PatchDataLayer.hpp:35
shamrock::solvergraph::IFieldRefs::get_field
PatchDataField< T > & get_field(u64 id) const
Get the underlying PatchDataField at the given id.
Definition IFieldRefs.hpp:54
forces.hpp
GSPH force computation using Riemann solver results.
UpdateDerivs.hpp
GSPH derivative update module.
iterative.hpp
Iterative Riemann solver for GSPH (van Leer 1997)
sham::inv_sat_positive
T inv_sat_positive(T v, T minvsat=T{1e-9}, T satval=T{0.}) noexcept
inverse saturated (positive numbers only)
Definition math.hpp:841
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
shammath
namespace for math utility
Definition AABB.hpp:26
shamrock
namespace for the main framework
Definition __init__.py:1
sphkernels.hpp
sph kernels
shamrock::patch::Patch
Patch object that contain generic patch information.
Definition Patch.hpp:33
shamrock::patch::Patch::id_patch
u64 id_patch
unique key that identify the patch
Definition Patch.hpp:86
Generated on Wed May 13 2026 13:21:49 for Shamrock by doxygen 1.9.8