DDC 0.2.0
Loading...
Searching...
No Matches
spline_evaluator_2d.hpp
1// Copyright (C) The DDC development team, see COPYRIGHT.md file
2//
3// SPDX-License-Identifier: MIT
4
5#pragma once
6
7#include <array>
8#include <cstddef>
9#include <type_traits>
10
11#include <ddc/ddc.hpp>
12
13#include <Kokkos_Core.hpp>
14
15#include "integrals.hpp"
16#include "periodic_extrapolation_rule.hpp"
17
18namespace ddc {
19
20/**
21 * @brief A class to evaluate, differentiate or integrate a 2D spline function.
22 *
23 * A class which contains an operator () which can be used to evaluate, differentiate or integrate a 2D spline function.
24 *
25 * @tparam ExecSpace The Kokkos execution space on which the spline evaluation is performed.
26 * @tparam MemorySpace The Kokkos memory space on which the data (spline coefficients and evaluation) is stored.
27 * @tparam BSplines1 The discrete dimension representing the B-splines along the first dimension of interest.
28 * @tparam BSplines2 The discrete dimension representing the B-splines along the second dimension of interest.
29 * @tparam EvaluationDDim1 The first discrete dimension on which evaluation points are defined.
30 * @tparam EvaluationDDim2 The second discrete dimension on which evaluation points are defined.
31 * @tparam LowerExtrapolationRule1 The lower extrapolation rule type along first dimension of interest.
32 * @tparam UpperExtrapolationRule1 The upper extrapolation rule type along first dimension of interest.
33 * @tparam LowerExtrapolationRule2 The lower extrapolation rule type along second dimension of interest.
34 * @tparam UpperExtrapolationRule2 The upper extrapolation rule type along second dimension of interest.
35 * @tparam IDimX A variadic template of all the discrete dimensions forming the full space (EvaluationDDim1 + EvaluationDDim2 + batched dimensions).
36 */
37template <
38 class ExecSpace,
39 class MemorySpace,
40 class BSplines1,
41 class BSplines2,
42 class EvaluationDDim1,
43 class EvaluationDDim2,
44 class LowerExtrapolationRule1,
45 class UpperExtrapolationRule1,
46 class LowerExtrapolationRule2,
47 class UpperExtrapolationRule2,
48 class... IDimX>
50{
51private:
52 /**
53 * @brief Tag to indicate that the value of the spline should be evaluated.
54 */
55 struct eval_type
56 {
57 };
58
59 /**
60 * @brief Tag to indicate that derivative of the spline should be evaluated.
61 */
62 struct eval_deriv_type
63 {
64 };
65
66public:
67 /// @brief The type of the first evaluation continuous dimension used by this class.
68 using continuous_dimension_type1 = typename BSplines1::continuous_dimension_type;
69
70 /// @brief The type of the second evaluation continuous dimension used by this class.
71 using continuous_dimension_type2 = typename BSplines2::continuous_dimension_type;
72
73 /// @brief The type of the Kokkos execution space used by this class.
74 using exec_space = ExecSpace;
75
76 /// @brief The type of the Kokkos memory space used by this class.
77 using memory_space = MemorySpace;
78
79 /// @brief The type of the first discrete dimension of interest used by this class.
80 using evaluation_discrete_dimension_type1 = EvaluationDDim1;
81
82 /// @brief The type of the second discrete dimension of interest used by this class.
83 using evaluation_discrete_dimension_type2 = EvaluationDDim2;
84
85 /// @brief The discrete dimension representing the B-splines along first dimension.
86 using bsplines_type1 = BSplines1;
87
88 /// @brief The discrete dimension representing the B-splines along second dimension.
89 using bsplines_type2 = BSplines2;
90
91 /// @brief The type of the domain for the 1D evaluation mesh along first dimension used by this class.
92 using evaluation_domain_type1 = ddc::DiscreteDomain<evaluation_discrete_dimension_type1>;
93
94 /// @brief The type of the domain for the 1D evaluation mesh along second dimension used by this class.
95 using evaluation_domain_type2 = ddc::DiscreteDomain<evaluation_discrete_dimension_type2>;
96
97 /// @brief The type of the domain for the 2D evaluation mesh used by this class.
98 using evaluation_domain_type = ddc::DiscreteDomain<
99 evaluation_discrete_dimension_type1,
100 evaluation_discrete_dimension_type2>;
101
102 /// @brief The type of the whole domain representing evaluation points.
103 using batched_evaluation_domain_type = ddc::DiscreteDomain<IDimX...>;
104
105 /// @brief The type of the 1D spline domain corresponding to the first dimension of interest.
106 using spline_domain_type1 = ddc::DiscreteDomain<bsplines_type1>;
107
108 /// @brief The type of the 1D spline domain corresponding to the second dimension of interest.
109 using spline_domain_type2 = ddc::DiscreteDomain<bsplines_type2>;
110
111 /// @brief The type of the 2D spline domain corresponding to the dimensions of interest.
112 using spline_domain_type = ddc::DiscreteDomain<bsplines_type1, bsplines_type2>;
113
114 /**
115 * @brief The type of the batch domain (obtained by removing the dimensions of interest
116 * from the whole domain).
117 */
118 using batch_domain_type = typename ddc::remove_dims_of_t<
119 batched_evaluation_domain_type,
120 evaluation_discrete_dimension_type1,
121 evaluation_discrete_dimension_type2>;
122
123 /**
124 * @brief The type of the whole spline domain (cartesian product of 2D spline domain
125 * and batch domain) preserving the underlying memory layout (order of dimensions).
126 */
127 using batched_spline_domain_type =
128 typename ddc::detail::convert_type_seq_to_discrete_domain_t<ddc::type_seq_replace_t<
129 ddc::detail::TypeSeq<IDimX...>,
130 ddc::detail::TypeSeq<
131 evaluation_discrete_dimension_type1,
132 evaluation_discrete_dimension_type2>,
133 ddc::detail::TypeSeq<bsplines_type1, bsplines_type2>>>;
134
135 /// @brief The type of the extrapolation rule at the lower boundary along the first dimension.
136 using lower_extrapolation_rule_1_type = LowerExtrapolationRule1;
137
138 /// @brief The type of the extrapolation rule at the upper boundary along the first dimension.
139 using upper_extrapolation_rule_1_type = UpperExtrapolationRule1;
140
141 /// @brief The type of the extrapolation rule at the lower boundary along the second dimension.
142 using lower_extrapolation_rule_2_type = LowerExtrapolationRule2;
143
144 /// @brief The type of the extrapolation rule at the upper boundary along the second dimension.
145 using upper_extrapolation_rule_2_type = UpperExtrapolationRule2;
146
147private:
148 LowerExtrapolationRule1 m_lower_extrap_rule_1;
149
150 UpperExtrapolationRule1 m_upper_extrap_rule_1;
151
152 LowerExtrapolationRule2 m_lower_extrap_rule_2;
153
154 UpperExtrapolationRule2 m_upper_extrap_rule_2;
155
156public:
157 static_assert(
158 std::is_same_v<LowerExtrapolationRule1,
159 typename ddc::PeriodicExtrapolationRule<continuous_dimension_type1>>
160 == bsplines_type1::is_periodic()
161 && std::is_same_v<
162 UpperExtrapolationRule1,
163 typename ddc::PeriodicExtrapolationRule<continuous_dimension_type1>>
164 == bsplines_type1::is_periodic()
165 && std::is_same_v<
166 LowerExtrapolationRule2,
167 typename ddc::PeriodicExtrapolationRule<continuous_dimension_type2>>
168 == bsplines_type2::is_periodic()
169 && std::is_same_v<
170 UpperExtrapolationRule2,
171 typename ddc::PeriodicExtrapolationRule<continuous_dimension_type2>>
172 == bsplines_type2::is_periodic(),
173 "PeriodicExtrapolationRule has to be used if and only if dimension is periodic");
174 static_assert(
175 std::is_invocable_r_v<
176 double,
177 LowerExtrapolationRule1,
178 ddc::Coordinate<continuous_dimension_type1>,
179 ddc::ChunkSpan<
180 double const,
181 spline_domain_type,
182 Kokkos::layout_right,
183 memory_space>>,
184 "LowerExtrapolationRule1::operator() has to be callable "
185 "with usual arguments.");
186 static_assert(
187 std::is_invocable_r_v<
188 double,
189 UpperExtrapolationRule1,
190 ddc::Coordinate<continuous_dimension_type1>,
191 ddc::ChunkSpan<
192 double const,
193 spline_domain_type,
194 Kokkos::layout_right,
195 memory_space>>,
196 "UpperExtrapolationRule1::operator() has to be callable "
197 "with usual arguments.");
198 static_assert(
199 std::is_invocable_r_v<
200 double,
201 LowerExtrapolationRule2,
202 ddc::Coordinate<continuous_dimension_type2>,
203 ddc::ChunkSpan<
204 double const,
205 spline_domain_type,
206 Kokkos::layout_right,
207 memory_space>>,
208 "LowerExtrapolationRule2::operator() has to be callable "
209 "with usual arguments.");
210 static_assert(
211 std::is_invocable_r_v<
212 double,
213 UpperExtrapolationRule2,
214 ddc::Coordinate<continuous_dimension_type2>,
215 ddc::ChunkSpan<
216 double const,
217 spline_domain_type,
218 Kokkos::layout_right,
219 memory_space>>,
220 "UpperExtrapolationRule2::operator() has to be callable "
221 "with usual arguments.");
222
223 /**
224 * @brief Build a SplineEvaluator2D acting on batched_spline_domain.
225 *
226 * @param lower_extrap_rule1 The extrapolation rule at the lower boundary along the first dimension.
227 * @param upper_extrap_rule1 The extrapolation rule at the upper boundary along the first dimension.
228 * @param lower_extrap_rule2 The extrapolation rule at the lower boundary along the second dimension.
229 * @param upper_extrap_rule2 The extrapolation rule at the upper boundary along the second dimension.
230 *
231 * @see NullExtrapolationRule ConstantExtrapolationRule PeriodicExtrapolationRule
232 */
233 explicit SplineEvaluator2D(
234 LowerExtrapolationRule1 const& lower_extrap_rule1,
235 UpperExtrapolationRule1 const& upper_extrap_rule1,
236 LowerExtrapolationRule2 const& lower_extrap_rule2,
237 UpperExtrapolationRule2 const& upper_extrap_rule2)
238 : m_lower_extrap_rule_1(lower_extrap_rule1)
239 , m_upper_extrap_rule_1(upper_extrap_rule1)
240 , m_lower_extrap_rule_2(lower_extrap_rule2)
241 , m_upper_extrap_rule_2(upper_extrap_rule2)
242 {
243 }
244
245 /**
246 * @brief Copy-constructs.
247 *
248 * @param x A reference to another SplineEvaluator.
249 */
250 SplineEvaluator2D(SplineEvaluator2D const& x) = default;
251
252 /**
253 * @brief Move-constructs.
254 *
255 * @param x An rvalue to another SplineEvaluator.
256 */
258
259 /// @brief Destructs.
260 ~SplineEvaluator2D() = default;
261
262 /**
263 * @brief Copy-assigns.
264 *
265 * @param x A reference to another SplineEvaluator.
266 * @return A reference to this object.
267 */
268 SplineEvaluator2D& operator=(SplineEvaluator2D const& x) = default;
269
270 /**
271 * @brief Move-assigns.
272 *
273 * @param x An rvalue to another SplineEvaluator.
274 * @return A reference to this object.
275 */
277
278 /**
279 * @brief Get the lower extrapolation rule along the first dimension.
280 *
281 * Extrapolation rules are functors used to define the behavior of the SplineEvaluator out of the domain where the break points of the B-splines are defined.
282 *
283 * @return The lower extrapolation rule along the first dimension.
284 *
285 * @see NullExtrapolationRule ConstantExtrapolationRule PeriodicExtrapolationRule
286 */
287 lower_extrapolation_rule_1_type lower_extrapolation_rule_dim_1() const
288 {
289 return m_lower_extrap_rule_1;
290 }
291
292 /**
293 * @brief Get the upper extrapolation rule along the first dimension.
294 *
295 * Extrapolation rules are functors used to define the behavior of the SplineEvaluator out of the domain where the break points of the B-splines are defined.
296 *
297 * @return The upper extrapolation rule along the first dimension.
298 *
299 * @see NullExtrapolationRule ConstantExtrapolationRule PeriodicExtrapolationRule
300 */
301 upper_extrapolation_rule_1_type upper_extrapolation_rule_dim_1() const
302 {
303 return m_upper_extrap_rule_1;
304 }
305
306 /**
307 * @brief Get the lower extrapolation rule along the second dimension.
308 *
309 * Extrapolation rules are functors used to define the behavior of the SplineEvaluator out of the domain where the break points of the B-splines are defined.
310 *
311 * @return The lower extrapolation rule along the second dimension.
312 *
313 * @see NullExtrapolationRule ConstantExtrapolationRule PeriodicExtrapolationRule
314 */
315 lower_extrapolation_rule_2_type lower_extrapolation_rule_dim_2() const
316 {
317 return m_lower_extrap_rule_2;
318 }
319
320 /**
321 * @brief Get the upper extrapolation rule along the second dimension.
322 *
323 * Extrapolation rules are functors used to define the behavior of the SplineEvaluator out of the domain where the break points of the B-splines are defined.
324 *
325 * @return The upper extrapolation rule along the second dimension.
326 *
327 * @see NullExtrapolationRule ConstantExtrapolationRule PeriodicExtrapolationRule
328 */
329 upper_extrapolation_rule_2_type upper_extrapolation_rule_dim_2() const
330 {
331 return m_upper_extrap_rule_2;
332 }
333
334 /**
335 * @brief Evaluate 2D spline function (described by its spline coefficients) at a given coordinate.
336 *
337 * The spline coefficients represent a 2D spline function defined on a B-splines (basis splines). They can be obtained via various methods, such as using a SplineBuilder2D.
338 *
339 * Remark: calling SplineBuilder2D then SplineEvaluator2D corresponds to a 2D spline interpolation.
340 *
341 * @param coord_eval The coordinate where the spline is evaluated. Note that only the components along the dimensions of interest are used.
342 * @param spline_coef A ChunkSpan storing the 2D spline coefficients.
343 *
344 * @return The value of the spline function at the desired coordinate.
345 */
346 template <class Layout, class... CoordsDims>
347 KOKKOS_FUNCTION double operator()(
348 ddc::Coordinate<CoordsDims...> const& coord_eval,
349 ddc::ChunkSpan<double const, spline_domain_type, Layout, memory_space> const
350 spline_coef) const
351 {
352 return eval(coord_eval, spline_coef);
353 }
354
355 /**
356 * @brief Evaluate 2D spline function (described by its spline coefficients) on a mesh.
357 *
358 * The spline coefficients represent a 2D spline function defined on a cartesian product of batch_domain and B-splines
359 * (basis splines). They can be obtained via various methods, such as using a SplineBuilder2D.
360 *
361 * This is not a nD evaluation. This is a batched 2D evaluation. This means that for each slice of coordinates
362 * identified by a batch_domain_type::discrete_element_type, the evaluation is performed with the 2D set of
363 * spline coefficients identified by the same batch_domain_type::discrete_element_type.
364 *
365 * Remark: calling SplineBuilder2D then SplineEvaluator2D corresponds to a 2D spline interpolation.
366 *
367 * @param[out] spline_eval The values of the 2D spline function at the desired coordinates. For practical reasons those are
368 * stored in a ChunkSpan defined on a batched_evaluation_domain_type.
369 * @param[in] coords_eval The coordinates where the spline is evaluated. Those are
370 * stored in a ChunkSpan defined on a batched_evaluation_domain_type. Note that the coordinates of the
371 * points represented by this domain are unused and irrelevant (but the points themselves (DiscreteElement) are used to select
372 * the set of 2D spline coefficients retained to perform the evaluation).
373 * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients.
374 */
375 template <class Layout1, class Layout2, class Layout3, class... CoordsDims>
376 void operator()(
377 ddc::ChunkSpan<double, batched_evaluation_domain_type, Layout1, memory_space> const
378 spline_eval,
379 ddc::ChunkSpan<
380 ddc::Coordinate<CoordsDims...> const,
381 batched_evaluation_domain_type,
382 Layout2,
383 memory_space> const coords_eval,
384 ddc::ChunkSpan<double const, batched_spline_domain_type, Layout3, memory_space> const
385 spline_coef) const
386 {
387 batch_domain_type const batch_domain(coords_eval.domain());
388 evaluation_domain_type1 const evaluation_domain1(spline_eval.domain());
389 evaluation_domain_type2 const evaluation_domain2(spline_eval.domain());
390 ddc::parallel_for_each(
391 "ddc_splines_evaluate_2d",
392 exec_space(),
393 batch_domain,
394 KOKKOS_CLASS_LAMBDA(typename batch_domain_type::discrete_element_type const j) {
395 const auto spline_eval_2D = spline_eval[j];
396 const auto coords_eval_2D = coords_eval[j];
397 const auto spline_coef_2D = spline_coef[j];
398 for (auto const i1 : evaluation_domain1) {
399 for (auto const i2 : evaluation_domain2) {
400 spline_eval_2D(i1, i2) = eval(coords_eval_2D(i1, i2), spline_coef_2D);
401 }
402 }
403 });
404 }
405
406 /**
407 * @brief Evaluate 2D spline function (described by its spline coefficients) on a mesh.
408 *
409 * The spline coefficients represent a 2D spline function defined on a cartesian product of batch_domain and B-splines
410 * (basis splines). They can be obtained via various methods, such as using a SplineBuilder2D.
411 *
412 * This is not a multidimensional evaluation. This is a batched 2D evaluation.
413 * This means that for each slice of spline_eval the evaluation is performed with
414 * the 2D set of spline coefficients identified by the same batch_domain_type::discrete_element_type.
415 *
416 * Remark: calling SplineBuilder2D then SplineEvaluator2D corresponds to a 2D spline interpolation.
417 *
418 * @param[out] spline_eval The values of the 2D spline function at their coordinates.
419 * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients.
420 */
421 template <class Layout1, class Layout2>
422 void operator()(
423 ddc::ChunkSpan<double, batched_evaluation_domain_type, Layout1, memory_space> const
424 spline_eval,
425 ddc::ChunkSpan<double const, batched_spline_domain_type, Layout2, memory_space> const
426 spline_coef) const
427 {
428 batch_domain_type const batch_domain(spline_eval.domain());
429 evaluation_domain_type1 const evaluation_domain1(spline_eval.domain());
430 evaluation_domain_type2 const evaluation_domain2(spline_eval.domain());
431 ddc::parallel_for_each(
432 "ddc_splines_evaluate_2d",
433 exec_space(),
434 batch_domain,
435 KOKKOS_CLASS_LAMBDA(typename batch_domain_type::discrete_element_type const j) {
436 const auto spline_eval_2D = spline_eval[j];
437 const auto spline_coef_2D = spline_coef[j];
438 for (auto const i1 : evaluation_domain1) {
439 for (auto const i2 : evaluation_domain2) {
440 ddc::Coordinate<continuous_dimension_type1, continuous_dimension_type2>
441 coord_eval_2D(ddc::coordinate(i1), ddc::coordinate(i2));
442 spline_eval_2D(i1, i2) = eval(coord_eval_2D(i1, i2), spline_coef_2D);
443 }
444 }
445 });
446 }
447
448 /**
449 * @brief Differentiate 2D spline function (described by its spline coefficients) at a given coordinate along first dimension of interest.
450 *
451 * The spline coefficients represent a 2D spline function defined on a B-splines (basis splines). They can be
452 * obtained via various methods, such as using a SplineBuilder2D.
453 *
454 * @param coord_eval The coordinate where the spline is differentiated. Note that only the components along the dimensions of interest are used.
455 * @param spline_coef A ChunkSpan storing the 2D spline coefficients.
456 *
457 * @return The derivative of the spline function at the desired coordinate.
458 */
459 template <class Layout, class... CoordsDims>
461 ddc::Coordinate<CoordsDims...> const& coord_eval,
462 ddc::ChunkSpan<double const, spline_domain_type, Layout, memory_space> const
463 spline_coef) const
464 {
465 return eval_no_bc<eval_deriv_type, eval_type>(coord_eval, spline_coef);
466 }
467
468 /**
469 * @brief Differentiate 2D spline function (described by its spline coefficients) at a given coordinate along second dimension of interest.
470 *
471 * The spline coefficients represent a 2D spline function defined on a B-splines (basis splines). They can be
472 * obtained via various methods, such as using a SplineBuilder2D.
473 *
474 * @param coord_eval The coordinate where the spline is differentiated. Note that only the components along the dimensions of interest are used.
475 * @param spline_coef A ChunkSpan storing the 2D spline coefficients.
476 *
477 * @return The derivative of the spline function at the desired coordinate.
478 */
479 template <class Layout, class... CoordsDims>
481 ddc::Coordinate<CoordsDims...> const& coord_eval,
482 ddc::ChunkSpan<double const, spline_domain_type, Layout, memory_space> const
483 spline_coef) const
484 {
485 return eval_no_bc<eval_type, eval_deriv_type>(coord_eval, spline_coef);
486 }
487
488 /**
489 * @brief Cross-differentiate 2D spline function (described by its spline coefficients) at a given coordinate.
490 *
491 * The spline coefficients represent a 2D spline function defined on a B-splines (basis splines). They can be
492 * obtained via various methods, such as using a SplineBuilder2D.
493 *
494 * @param coord_eval The coordinate where the spline is differentiated. Note that only the components along the dimensions of interest are used.
495 * @param spline_coef A ChunkSpan storing the 2D spline coefficients.
496 *
497 * @return The derivative of the spline function at the desired coordinate.
498 */
499 template <class Layout, class... CoordsDims>
501 ddc::Coordinate<CoordsDims...> const& coord_eval,
502 ddc::ChunkSpan<double const, spline_domain_type, Layout, memory_space> const
503 spline_coef) const
504 {
505 return eval_no_bc<eval_deriv_type, eval_deriv_type>(coord_eval, spline_coef);
506 }
507
508 /**
509 * @brief Differentiate 2D spline function (described by its spline coefficients) at a given coordinate along a specified dimension of interest.
510 *
511 * The spline coefficients represent a 2D spline function defined on a B-splines (basis splines). They can be
512 * obtained via various methods, such as using a SplineBuilder2D.
513 *
514 * @tparam InterestDim Dimension along which differentiation is performed.
515 *
516 * @param coord_eval The coordinate where the spline is differentiated. Note that only the components along the dimensions of interest are used.
517 * @param spline_coef A ChunkSpan storing the 2D spline coefficients.
518 *
519 * @return The derivative of the spline function at the desired coordinate.
520 */
521 template <class InterestDim, class Layout, class... CoordsDims>
522 KOKKOS_FUNCTION double deriv(
523 ddc::Coordinate<CoordsDims...> const& coord_eval,
524 ddc::ChunkSpan<double const, spline_domain_type, Layout, memory_space> const
525 spline_coef) const
526 {
527 static_assert(
528 std::is_same_v<InterestDim, continuous_dimension_type1>
529 || std::is_same_v<InterestDim, continuous_dimension_type2>);
530 if constexpr (std::is_same_v<
531 InterestDim,
532 typename evaluation_discrete_dimension_type1::
533 continuous_dimension_type>) {
534 return deriv_dim_1(coord_eval, spline_coef);
535 } else if constexpr (std::is_same_v<
536 InterestDim,
537 typename evaluation_discrete_dimension_type2::
538 continuous_dimension_type>) {
539 return deriv_dim_2(coord_eval, spline_coef);
540 }
541 }
542
543 /**
544 * @brief Double-differentiate 2D spline function (described by its spline coefficients) at a given coordinate along specified dimensions of interest.
545 *
546 * The spline coefficients represent a 2D spline function defined on a B-splines (basis splines). They can be
547 * obtained via various methods, such as using a SplineBuilder2D.
548 *
549 * Note: double-differentiation other than cross-differentiation is not supported atm. See #440
550 *
551 * @tparam InterestDim1 First dimension along which differentiation is performed.
552 * @tparam InterestDim2 Second dimension along which differentiation is performed.
553 *
554 * @param coord_eval The coordinate where the spline is double-differentiated. Note that only the components along the dimensions of interest are used.
555 * @param spline_coef A ChunkSpan storing the 2D spline coefficients.
556 *
557 * @return The derivative of the spline function at the desired coordinate.
558 */
559 template <class InterestDim1, class InterestDim2, class Layout, class... CoordsDims>
560 KOKKOS_FUNCTION double deriv2(
561 ddc::Coordinate<CoordsDims...> const& coord_eval,
562 ddc::ChunkSpan<double const, spline_domain_type, Layout, memory_space> const
563 spline_coef) const
564 {
565 static_assert(
566 (std::is_same_v<
567 InterestDim1,
568 typename evaluation_discrete_dimension_type1::continuous_dimension_type>
569 && std::is_same_v<InterestDim2, continuous_dimension_type2>)
570 || (std::is_same_v<
571 InterestDim2,
572 typename evaluation_discrete_dimension_type1::continuous_dimension_type>
573 && std::is_same_v<InterestDim1, continuous_dimension_type2>));
574 return deriv_1_and_2(coord_eval, spline_coef);
575 }
576
577 /**
578 * @brief Differentiate 2D spline function (described by its spline coefficients) on a mesh along first dimension of interest.
579 *
580 * The spline coefficients represent a 2D spline function defined on a cartesian product of batch_domain and B-splines
581 * (basis splines). They can be obtained via various methods, such as using a SplineBuilder2D.
582 *
583 * This is not a nD evaluation. This is a batched 2D differentiation.
584 * This means that for each slice of coordinates identified by a batch_domain_type::discrete_element_type,
585 * the differentiation is performed with the 2D set of spline coefficients identified by the same batch_domain_type::discrete_element_type.
586 *
587 * @param[out] spline_eval The derivatives of the 2D spline function at the desired coordinates. For practical reasons those are
588 * stored in a ChunkSpan defined on a batched_evaluation_domain_type. Note that the coordinates of the
589 * points represented by this domain are unused and irrelevant (but the points themselves (DiscreteElement) are used to select
590 * the set of 2D spline coefficients retained to perform the evaluation).
591 * @param[in] coords_eval The coordinates where the spline is differentiated. Those are
592 * stored in a ChunkSpan defined on a batched_evaluation_domain_type. Note that the coordinates of the
593 * points represented by this domain are unused and irrelevant (but the points themselves (DiscreteElement) are used to select
594 * the set of 2D spline coefficients retained to perform the evaluation).
595 * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients.
596 */
597 template <class Layout1, class Layout2, class Layout3, class... CoordsDims>
598 void deriv_dim_1(
599 ddc::ChunkSpan<double, batched_evaluation_domain_type, Layout1, memory_space> const
600 spline_eval,
601 ddc::ChunkSpan<
602 ddc::Coordinate<CoordsDims...> const,
603 batched_evaluation_domain_type,
604 Layout2,
605 memory_space> const coords_eval,
606 ddc::ChunkSpan<double const, batched_spline_domain_type, Layout3, memory_space> const
607 spline_coef) const
608 {
609 batch_domain_type const batch_domain(coords_eval.domain());
610 evaluation_domain_type1 const evaluation_domain1(spline_eval.domain());
611 evaluation_domain_type2 const evaluation_domain2(spline_eval.domain());
612 ddc::parallel_for_each(
613 "ddc_splines_differentiate_2d_dim_1",
614 exec_space(),
615 batch_domain,
616 KOKKOS_CLASS_LAMBDA(typename batch_domain_type::discrete_element_type const j) {
617 const auto spline_eval_2D = spline_eval[j];
618 const auto coords_eval_2D = coords_eval[j];
619 const auto spline_coef_2D = spline_coef[j];
620 for (auto const i1 : evaluation_domain1) {
621 for (auto const i2 : evaluation_domain2) {
622 spline_eval_2D(i1, i2) = eval_no_bc<
623 eval_deriv_type,
624 eval_type>(coords_eval_2D(i1, i2), spline_coef_2D);
625 }
626 }
627 });
628 }
629
630 /**
631 * @brief Differentiate 2D spline function (described by its spline coefficients) on a mesh along first dimension of interest.
632 *
633 * The spline coefficients represent a 2D spline function defined on a cartesian product of batch_domain and B-splines
634 * (basis splines). They can be obtained via various methods, such as using a SplineBuilder2D.
635 *
636 * This is not a multidimensional evaluation. This is a batched 2D evaluation.
637 * This means that for each slice of spline_eval the evaluation is performed with
638 * the 2D set of spline coefficients identified by the same batch_domain_type::discrete_element_type.
639 *
640 * @param[out] spline_eval The derivatives of the 2D spline function at the desired coordinates.
641 * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients.
642 */
643 template <class Layout1, class Layout2>
644 void deriv_dim_1(
645 ddc::ChunkSpan<double, batched_evaluation_domain_type, Layout1, memory_space> const
646 spline_eval,
647 ddc::ChunkSpan<double const, batched_spline_domain_type, Layout2, memory_space> const
648 spline_coef) const
649 {
650 batch_domain_type const batch_domain(spline_eval.domain());
651 evaluation_domain_type1 const evaluation_domain1(spline_eval.domain());
652 evaluation_domain_type2 const evaluation_domain2(spline_eval.domain());
653 ddc::parallel_for_each(
654 "ddc_splines_differentiate_2d_dim_1",
655 exec_space(),
656 batch_domain,
657 KOKKOS_CLASS_LAMBDA(typename batch_domain_type::discrete_element_type const j) {
658 const auto spline_eval_2D = spline_eval[j];
659 const auto spline_coef_2D = spline_coef[j];
660 for (auto const i1 : evaluation_domain1) {
661 for (auto const i2 : evaluation_domain2) {
662 ddc::Coordinate<continuous_dimension_type1, continuous_dimension_type2>
663 coord_eval_2D(ddc::coordinate(i1), ddc::coordinate(i2));
664 spline_eval_2D(i1, i2) = eval_no_bc<
665 eval_deriv_type,
666 eval_type>(coord_eval_2D, spline_coef_2D);
667 }
668 }
669 });
670 }
671
672 /**
673 * @brief Differentiate 2D spline function (described by its spline coefficients) on a mesh along second dimension of interest.
674 *
675 * The spline coefficients represent a 2D spline function defined on a cartesian product of batch_domain and B-splines
676 * (basis splines). They can be obtained via various methods, such as using a SplineBuilder2D.
677 *
678 * This is not a nD differentiation. This is a batched 2D differentiation.
679 * This means that for each slice of coordinates identified by a batch_domain_type::discrete_element_type,
680 * the differentiation is performed with the 2D set of spline coefficients identified by the same batch_domain_type::discrete_element_type.
681 *
682 * @param[out] spline_eval The derivatives of the 2D spline function at the desired coordinates. For practical reasons those are
683 * stored in a ChunkSpan defined on a batched_evaluation_domain_type.
684 * @param[in] coords_eval The coordinates where the spline is differentiated. Those are
685 * stored in a ChunkSpan defined on a batched_evaluation_domain_type. Note that the coordinates of the
686 * points represented by this domain are unused and irrelevant (but the points themselves (DiscreteElement) are used to select
687 * the set of 2D spline coefficients retained to perform the evaluation).
688 * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients.
689 */
690 template <class Layout1, class Layout2, class Layout3, class... CoordsDims>
691 void deriv_dim_2(
692 ddc::ChunkSpan<double, batched_evaluation_domain_type, Layout1, memory_space> const
693 spline_eval,
694 ddc::ChunkSpan<
695 ddc::Coordinate<CoordsDims...> const,
696 batched_evaluation_domain_type,
697 Layout2,
698 memory_space> const coords_eval,
699 ddc::ChunkSpan<double const, batched_spline_domain_type, Layout3, memory_space> const
700 spline_coef) const
701 {
702 batch_domain_type const batch_domain(coords_eval.domain());
703 evaluation_domain_type1 const evaluation_domain1(spline_eval.domain());
704 evaluation_domain_type2 const evaluation_domain2(spline_eval.domain());
705 ddc::parallel_for_each(
706 "ddc_splines_differentiate_2d_dim_2",
707 exec_space(),
708 batch_domain,
709 KOKKOS_CLASS_LAMBDA(typename batch_domain_type::discrete_element_type const j) {
710 const auto spline_eval_2D = spline_eval[j];
711 const auto coords_eval_2D = coords_eval[j];
712 const auto spline_coef_2D = spline_coef[j];
713 for (auto const i1 : evaluation_domain1) {
714 for (auto const i2 : evaluation_domain2) {
715 spline_eval_2D(i1, i2) = eval_no_bc<
716 eval_type,
717 eval_deriv_type>(coords_eval_2D(i1, i2), spline_coef_2D);
718 }
719 }
720 });
721 }
722
723 /**
724 * @brief Differentiate 2D spline function (described by its spline coefficients) on a mesh along second dimension of interest.
725 *
726 * The spline coefficients represent a 2D spline function defined on a cartesian product of batch_domain and B-splines
727 * (basis splines). They can be obtained via various methods, such as using a SplineBuilder2D.
728 *
729 * This is not a multidimensional evaluation. This is a batched 2D evaluation.
730 * This means that for each slice of spline_eval the evaluation is performed with
731 * the 2D set of spline coefficients identified by the same batch_domain_type::discrete_element_type.
732 *
733 * @param[out] spline_eval The derivatives of the 2D spline function at the desired coordinates.
734 * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients.
735 */
736 template <class Layout1, class Layout2>
737 void deriv_dim_2(
738 ddc::ChunkSpan<double, batched_evaluation_domain_type, Layout1, memory_space> const
739 spline_eval,
740 ddc::ChunkSpan<double const, batched_spline_domain_type, Layout2, memory_space> const
741 spline_coef) const
742 {
743 batch_domain_type const batch_domain(spline_eval.domain());
744 evaluation_domain_type1 const evaluation_domain1(spline_eval.domain());
745 evaluation_domain_type2 const evaluation_domain2(spline_eval.domain());
746 ddc::parallel_for_each(
747 "ddc_splines_differentiate_2d_dim_2",
748 exec_space(),
749 batch_domain,
750 KOKKOS_CLASS_LAMBDA(typename batch_domain_type::discrete_element_type const j) {
751 const auto spline_eval_2D = spline_eval[j];
752 const auto spline_coef_2D = spline_coef[j];
753 for (auto const i1 : evaluation_domain1) {
754 for (auto const i2 : evaluation_domain2) {
755 ddc::Coordinate<continuous_dimension_type1, continuous_dimension_type2>
756 coord_eval_2D(ddc::coordinate(i1), ddc::coordinate(i2));
757 spline_eval_2D(i1, i2) = eval_no_bc<
758 eval_type,
759 eval_deriv_type>(coord_eval_2D, spline_coef_2D);
760 }
761 }
762 });
763 }
764
765 /**
766 * @brief Cross-differentiate 2D spline function (described by its spline coefficients) on a mesh along dimensions of interest.
767 *
768 * The spline coefficients represent a 2D spline function defined on a cartesian product of batch_domain and B-splines
769 * (basis splines). They can be obtained via various methods, such as using a SplineBuilder2D.
770 *
771 * This is not a nD cross-differentiation. This is a batched 2D cross-differentiation.
772 * This means that for each slice of coordinates identified by a batch_domain_type::discrete_element_type,
773 * the cross-differentiation is performed with the 2D set of spline coefficients identified by the same batch_domain_type::discrete_element_type.
774 *
775 * @param[out] spline_eval The cross-derivatives of the 2D spline function at the desired coordinates. For practical reasons those are
776 * stored in a ChunkSpan defined on a batched_evaluation_domain_type.
777 * @param[in] coords_eval The coordinates where the spline is differentiated. Those are
778 * stored in a ChunkSpan defined on a batched_evaluation_domain_type. Note that the coordinates of the
779 * points represented by this domain are unused and irrelevant (but the points themselves (DiscreteElement) are used to select
780 * the set of 2D spline coefficients retained to perform the evaluation).
781 * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients.
782 */
783 template <class Layout1, class Layout2, class Layout3, class... CoordsDims>
784 void deriv_1_and_2(
785 ddc::ChunkSpan<double, batched_evaluation_domain_type, Layout1, memory_space> const
786 spline_eval,
787 ddc::ChunkSpan<
788 ddc::Coordinate<CoordsDims...> const,
789 batched_evaluation_domain_type,
790 Layout2,
791 memory_space> const coords_eval,
792 ddc::ChunkSpan<double const, batched_spline_domain_type, Layout3, memory_space> const
793 spline_coef) const
794 {
795 batch_domain_type const batch_domain(coords_eval.domain());
796 evaluation_domain_type1 const evaluation_domain1(spline_eval.domain());
797 evaluation_domain_type2 const evaluation_domain2(spline_eval.domain());
798 ddc::parallel_for_each(
799 "ddc_splines_cross_differentiate",
800 exec_space(),
801 batch_domain,
802 KOKKOS_CLASS_LAMBDA(typename batch_domain_type::discrete_element_type const j) {
803 const auto spline_eval_2D = spline_eval[j];
804 const auto coords_eval_2D = coords_eval[j];
805 const auto spline_coef_2D = spline_coef[j];
806 for (auto const i1 : evaluation_domain1) {
807 for (auto const i2 : evaluation_domain2) {
808 spline_eval_2D(i1, i2) = eval_no_bc<
809 eval_deriv_type,
810 eval_deriv_type>(coords_eval_2D(i1, i2), spline_coef_2D);
811 }
812 }
813 });
814 }
815
816 /**
817 * @brief Cross-differentiate 2D spline function (described by its spline coefficients) on a mesh along dimensions of interest.
818 *
819 * The spline coefficients represent a 2D spline function defined on a cartesian product of batch_domain and B-splines
820 * (basis splines). They can be obtained via various methods, such as using a SplineBuilder2D.
821 *
822 * This is not a multidimensional evaluation. This is a batched 2D evaluation.
823 * This means that for each slice of spline_eval the evaluation is performed with
824 * the 2D set of spline coefficients identified by the same batch_domain_type::discrete_element_type.
825 *
826 * @param[out] spline_eval The cross-derivatives of the 2D spline function at the desired coordinates.
827 * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients.
828 */
829 template <class Layout1, class Layout2>
830 void deriv_1_and_2(
831 ddc::ChunkSpan<double, batched_evaluation_domain_type, Layout1, memory_space> const
832 spline_eval,
833 ddc::ChunkSpan<double const, batched_spline_domain_type, Layout2, memory_space> const
834 spline_coef) const
835 {
836 batch_domain_type const batch_domain(spline_eval.domain());
837 evaluation_domain_type1 const evaluation_domain1(spline_eval.domain());
838 evaluation_domain_type2 const evaluation_domain2(spline_eval.domain());
839 ddc::parallel_for_each(
840 "ddc_splines_cross_differentiate",
841 exec_space(),
842 batch_domain,
843 KOKKOS_CLASS_LAMBDA(typename batch_domain_type::discrete_element_type const j) {
844 const auto spline_eval_2D = spline_eval[j];
845 const auto spline_coef_2D = spline_coef[j];
846 for (auto const i1 : evaluation_domain1) {
847 for (auto const i2 : evaluation_domain2) {
848 ddc::Coordinate<continuous_dimension_type1, continuous_dimension_type2>
849 coord_eval_2D(ddc::coordinate(i1), ddc::coordinate(i2));
850 spline_eval_2D(i1, i2) = eval_no_bc<
851 eval_deriv_type,
852 eval_deriv_type>(coord_eval_2D, spline_coef_2D);
853 }
854 }
855 });
856 }
857
858 /**
859 * @brief Differentiate spline function (described by its spline coefficients) on a mesh along a specified dimension of interest.
860 *
861 * The spline coefficients represent a 2D spline function defined on a cartesian product of batch_domain and B-splines
862 * (basis splines). They can be obtained via various methods, such as using a SplineBuilder2D.
863 *
864 * This is not a nD evaluation. This is a batched 2D differentiation.
865 * This means that for each slice of coordinates identified by a batch_domain_type::discrete_element_type,
866 * the differentiation is performed with the 2D set of spline coefficients identified by the same batch_domain_type::discrete_element_type.
867 *
868 * @tparam InterestDim Dimension along which differentiation is performed.
869 * @param[out] spline_eval The derivatives of the 2D spline function at the desired coordinates. For practical reasons those are
870 * stored in a ChunkSpan defined on a batched_evaluation_domain_type.
871 * @param[in] coords_eval The coordinates where the spline is differentiated. Those are
872 * stored in a ChunkSpan defined on a batched_evaluation_domain_type. Note that the coordinates of the
873 * points represented by this domain are unused and irrelevant (but the points themselves (DiscreteElement) are used to select
874 * the set of 2D spline coefficients retained to perform the evaluation).
875 * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients.
876 */
877 template <class InterestDim, class Layout1, class Layout2, class Layout3, class... CoordsDims>
878 void deriv(
879 ddc::ChunkSpan<double, batched_evaluation_domain_type, Layout1, memory_space> const
880 spline_eval,
881 ddc::ChunkSpan<
882 ddc::Coordinate<CoordsDims...> const,
883 batched_evaluation_domain_type,
884 Layout2,
885 memory_space> const coords_eval,
886 ddc::ChunkSpan<double const, batched_spline_domain_type, Layout3, memory_space> const
887 spline_coef) const
888 {
889 static_assert(
890 std::is_same_v<
891 InterestDim,
892 typename evaluation_discrete_dimension_type1::continuous_dimension_type>
893 || std::is_same_v<InterestDim, continuous_dimension_type2>);
894 if constexpr (std::is_same_v<
895 InterestDim,
896 typename evaluation_discrete_dimension_type1::
897 continuous_dimension_type>) {
898 return deriv_dim_1(spline_eval, coords_eval, spline_coef);
899 } else if constexpr (std::is_same_v<
900 InterestDim,
901 typename evaluation_discrete_dimension_type2::
902 continuous_dimension_type>) {
903 return deriv_dim_2(spline_eval, coords_eval, spline_coef);
904 }
905 }
906
907 /**
908 * @brief Differentiate spline function (described by its spline coefficients) on a mesh along a specified dimension of interest.
909 *
910 * The spline coefficients represent a 2D spline function defined on a cartesian product of batch_domain and B-splines
911 * (basis splines). They can be obtained via various methods, such as using a SplineBuilder2D.
912 *
913 * This is not a multidimensional evaluation. This is a batched 2D evaluation.
914 * This means that for each slice of spline_eval the evaluation is performed with
915 * the 2D set of spline coefficients identified by the same batch_domain_type::discrete_element_type.
916 *
917 * @tparam InterestDim Dimension along which differentiation is performed.
918 * @param[out] spline_eval The derivatives of the 2D spline function at the desired coordinates.
919 * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients.
920 */
921 template <class InterestDim, class Layout1, class Layout2>
922 void deriv(
923 ddc::ChunkSpan<double, batched_evaluation_domain_type, Layout1, memory_space> const
924 spline_eval,
925 ddc::ChunkSpan<double const, batched_spline_domain_type, Layout2, memory_space> const
926 spline_coef) const
927 {
928 static_assert(
929 std::is_same_v<
930 InterestDim,
931 typename evaluation_discrete_dimension_type1::continuous_dimension_type>
932 || std::is_same_v<InterestDim, continuous_dimension_type2>);
933 if constexpr (std::is_same_v<
934 InterestDim,
935 typename evaluation_discrete_dimension_type1::
936 continuous_dimension_type>) {
937 return deriv_dim_1(spline_eval, spline_coef);
938 } else if constexpr (std::is_same_v<
939 InterestDim,
940 typename evaluation_discrete_dimension_type2::
941 continuous_dimension_type>) {
942 return deriv_dim_2(spline_eval, spline_coef);
943 }
944 }
945
946 /**
947 * @brief Double-differentiate 2D spline function (described by its spline coefficients) on a mesh along specified dimensions of interest.
948 *
949 * The spline coefficients represent a 2D spline function defined on a cartesian product of batch_domain and B-splines
950 * (basis splines). They can be obtained via various methods, such as using a SplineBuilder2D.
951 *
952 * This is not a nD evaluation. This is a batched 2D differentiation.
953 * This means that for each slice of coordinates identified by a batch_domain_type::discrete_element_type,
954 * the differentiation is performed with the 2D set of spline coefficients identified by the same batch_domain_type::discrete_element_type.
955 *
956 * Note: double-differentiation other than cross-differentiation is not supported atm. See #440
957 *
958 * @tparam InterestDim1 First dimension along which differentiation is performed.
959 * @tparam InterestDim2 Second dimension along which differentiation is performed.
960 *
961 * @param[out] spline_eval The derivatives of the 2D spline function at the desired coordinates. For practical reasons those are
962 * stored in a ChunkSpan defined on a batched_evaluation_domain_type.
963 * @param[in] coords_eval The coordinates where the spline is differentiated. Those are
964 * stored in a ChunkSpan defined on a batched_evaluation_domain_type. Note that the coordinates of the
965 * points represented by this domain are unused and irrelevant (but the points themselves (DiscreteElement) are used to select
966 * the set of 2D spline coefficients retained to perform the evaluation).
967 * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients.
968 */
969 template <
970 class InterestDim1,
971 class InterestDim2,
972 class Layout1,
973 class Layout2,
974 class Layout3,
975 class... CoordsDims>
976 void deriv2(
977 ddc::ChunkSpan<double, batched_evaluation_domain_type, Layout1, memory_space> const
978 spline_eval,
979 ddc::ChunkSpan<
980 ddc::Coordinate<CoordsDims...> const,
981 batched_evaluation_domain_type,
982 Layout2,
983 memory_space> const coords_eval,
984 ddc::ChunkSpan<double const, batched_spline_domain_type, Layout3, memory_space> const
985 spline_coef) const
986 {
987 static_assert(
988 (std::is_same_v<
989 InterestDim1,
990 typename evaluation_discrete_dimension_type1::continuous_dimension_type>
991 && std::is_same_v<InterestDim2, continuous_dimension_type2>)
992 || (std::is_same_v<
993 InterestDim2,
994 typename evaluation_discrete_dimension_type1::continuous_dimension_type>
995 && std::is_same_v<InterestDim1, continuous_dimension_type2>));
996 return deriv_1_and_2(spline_eval, coords_eval, spline_coef);
997 }
998
999 /**
1000 * @brief Double-differentiate 2D spline function (described by its spline coefficients) on a mesh along specified dimensions of interest.
1001 *
1002 * The spline coefficients represent a 2D spline function defined on a cartesian product of batch_domain and B-splines
1003 * (basis splines). They can be obtained via various methods, such as using a SplineBuilder2D.
1004 *
1005 * This is not a multidimensional evaluation. This is a batched 2D evaluation.
1006 * This means that for each slice of spline_eval the evaluation is performed with
1007 * the 2D set of spline coefficients identified by the same batch_domain_type::discrete_element_type.
1008 *
1009 * Note: double-differentiation other than cross-differentiation is not supported atm. See #440
1010 *
1011 * @tparam InterestDim1 First dimension along which differentiation is performed.
1012 * @tparam InterestDim2 Second dimension along which differentiation is performed.
1013 *
1014 * @param[out] spline_eval The derivatives of the 2D spline function at the desired coordinates.
1015 * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients.
1016 */
1017 template <class InterestDim1, class InterestDim2, class Layout1, class Layout2>
1018 void deriv2(
1019 ddc::ChunkSpan<double, batched_evaluation_domain_type, Layout1, memory_space> const
1020 spline_eval,
1021 ddc::ChunkSpan<double const, batched_spline_domain_type, Layout2, memory_space> const
1022 spline_coef) const
1023 {
1024 static_assert(
1025 (std::is_same_v<
1026 InterestDim1,
1027 typename evaluation_discrete_dimension_type1::continuous_dimension_type>
1028 && std::is_same_v<InterestDim2, continuous_dimension_type2>)
1029 || (std::is_same_v<
1030 InterestDim2,
1031 typename evaluation_discrete_dimension_type1::continuous_dimension_type>
1032 && std::is_same_v<InterestDim1, continuous_dimension_type2>));
1033 return deriv_1_and_2(spline_eval, spline_coef);
1034 }
1035
1036 /** @brief Perform batched 2D integrations of a spline function (described by its spline coefficients) along the dimensions of interest and store results on a subdomain of batch_domain.
1037 *
1038 * The spline coefficients represent a 2D spline function defined on a B-splines (basis splines). They can be obtained via various methods, such as using a SplineBuilder2D.
1039 *
1040 * This is not a nD integration. This is a batched 2D integration.
1041 * This means that for each element of integrals, the integration is performed with the 2D set of
1042 * spline coefficients identified by the same DiscreteElement.
1043 *
1044 * @param[out] integrals The integrals of the 2D spline function on the subdomain of batch_domain. For practical reasons those are
1045 * stored in a ChunkSpan defined on a batch_domain_type. Note that the coordinates of the
1046 * points represented by this domain are unused and irrelevant.
1047 * @param[in] spline_coef A ChunkSpan storing the 2D spline coefficients.
1048 */
1049 template <class Layout1, class Layout2>
1050 void integrate(
1051 ddc::ChunkSpan<double, batch_domain_type, Layout1, memory_space> const integrals,
1052 ddc::ChunkSpan<double const, batched_spline_domain_type, Layout2, memory_space> const
1053 spline_coef) const
1054 {
1055 batch_domain_type batch_domain(integrals.domain());
1056 ddc::Chunk values1_alloc(
1057 ddc::DiscreteDomain<bsplines_type1>(spline_coef.domain()),
1058 ddc::KokkosAllocator<double, memory_space>());
1059 ddc::ChunkSpan values1 = values1_alloc.span_view();
1060 ddc::integrals(exec_space(), values1);
1061 ddc::Chunk values2_alloc(
1062 ddc::DiscreteDomain<bsplines_type2>(spline_coef.domain()),
1063 ddc::KokkosAllocator<double, memory_space>());
1064 ddc::ChunkSpan values2 = values2_alloc.span_view();
1065 ddc::integrals(exec_space(), values2);
1066
1067 ddc::parallel_for_each(
1068 "ddc_splines_integrate_bsplines",
1069 exec_space(),
1070 batch_domain,
1071 KOKKOS_LAMBDA(typename batch_domain_type::discrete_element_type const j) {
1072 integrals(j) = 0;
1073 for (typename spline_domain_type1::discrete_element_type const i1 :
1074 values1.domain()) {
1075 for (typename spline_domain_type2::discrete_element_type const i2 :
1076 values2.domain()) {
1077 integrals(j) += spline_coef(i1, i2, j) * values1(i1) * values2(i2);
1078 }
1079 }
1080 });
1081 }
1082
1083private:
1084 /**
1085 * @brief Evaluate the function on B-splines at the coordinate given.
1086 *
1087 * This function firstly deals with the boundary conditions and calls the SplineEvaluator2D::eval_no_bc function
1088 * to evaluate.
1089 *
1090 * @param[in] coord_eval
1091 * The 2D coordinate where we want to evaluate.
1092 * @param[in] spline_coef
1093 * The B-splines coefficients of the function we want to evaluate.
1094 * @param[out] vals1
1095 * A ChunkSpan with the not-null values of each function of the spline in the first dimension.
1096 * @param[out] vals2
1097 * A ChunkSpan with the not-null values of each function of the spline in the second dimension.
1098 *
1099 * @return A double with the value of the function at the coordinate given.
1100 *
1101 * @see SplineBoundaryValue
1102 */
1103 template <class Layout, class... CoordsDims>
1104 KOKKOS_INLINE_FUNCTION double eval(
1105 ddc::Coordinate<CoordsDims...> coord_eval,
1106 ddc::ChunkSpan<double const, spline_domain_type, Layout, memory_space> const
1107 spline_coef) const
1108 {
1109 using Dim1 = continuous_dimension_type1;
1110 using Dim2 = continuous_dimension_type2;
1111 if constexpr (bsplines_type1::is_periodic()) {
1112 if (ddc::get<Dim1>(coord_eval) < ddc::discrete_space<bsplines_type1>().rmin()
1113 || ddc::get<Dim1>(coord_eval) > ddc::discrete_space<bsplines_type1>().rmax()) {
1114 ddc::get<Dim1>(coord_eval)
1115 -= Kokkos::floor(
1116 (ddc::get<Dim1>(coord_eval)
1117 - ddc::discrete_space<bsplines_type1>().rmin())
1118 / ddc::discrete_space<bsplines_type1>().length())
1119 * ddc::discrete_space<bsplines_type1>().length();
1120 }
1121 }
1122 if constexpr (bsplines_type2::is_periodic()) {
1123 if (ddc::get<Dim2>(coord_eval) < ddc::discrete_space<bsplines_type2>().rmin()
1124 || ddc::get<Dim2>(coord_eval) > ddc::discrete_space<bsplines_type2>().rmax()) {
1125 ddc::get<Dim2>(coord_eval)
1126 -= Kokkos::floor(
1127 (ddc::get<Dim2>(coord_eval)
1128 - ddc::discrete_space<bsplines_type2>().rmin())
1129 / ddc::discrete_space<bsplines_type2>().length())
1130 * ddc::discrete_space<bsplines_type2>().length();
1131 }
1132 }
1133 if constexpr (!bsplines_type1::is_periodic()) {
1134 if (ddc::get<Dim1>(coord_eval) < ddc::discrete_space<bsplines_type1>().rmin()) {
1135 return m_lower_extrap_rule_1(coord_eval, spline_coef);
1136 }
1137 if (ddc::get<Dim1>(coord_eval) > ddc::discrete_space<bsplines_type1>().rmax()) {
1138 return m_upper_extrap_rule_1(coord_eval, spline_coef);
1139 }
1140 }
1141 if constexpr (!bsplines_type2::is_periodic()) {
1142 if (ddc::get<Dim2>(coord_eval) < ddc::discrete_space<bsplines_type2>().rmin()) {
1143 return m_lower_extrap_rule_2(coord_eval, spline_coef);
1144 }
1145 if (ddc::get<Dim2>(coord_eval) > ddc::discrete_space<bsplines_type2>().rmax()) {
1146 return m_upper_extrap_rule_2(coord_eval, spline_coef);
1147 }
1148 }
1149 return eval_no_bc<eval_type, eval_type>(
1150 ddc::Coordinate<continuous_dimension_type1, continuous_dimension_type2>(
1151 ddc::get<Dim1>(coord_eval),
1152 ddc::get<Dim2>(coord_eval)),
1153 spline_coef);
1154 }
1155
1156 /**
1157 * @brief Evaluate the function or its derivative at the coordinate given.
1158 *
1159 * @param[in] coord_eval
1160 * The coordinate where we want to evaluate.
1161 * @param[in] splne_coef
1162 * The B-splines coefficients of the function we want to evaluate.
1163 * @tparam EvalType1
1164 * A flag indicating if we evaluate the function or its derivative in the first dimension.
1165 * The type of this object is either `eval_type` or `eval_deriv_type`.
1166 * @tparam EvalType2
1167 * A flag indicating if we evaluate the function or its derivative in the second dimension.
1168 * The type of this object is either `eval_type` or `eval_deriv_type`.
1169 */
1170 template <class EvalType1, class EvalType2, class Layout, class... CoordsDims>
1171 KOKKOS_INLINE_FUNCTION double eval_no_bc(
1172 ddc::Coordinate<CoordsDims...> const& coord_eval,
1173 ddc::ChunkSpan<double const, spline_domain_type, Layout, memory_space> const
1174 spline_coef) const
1175 {
1176 static_assert(
1177 std::is_same_v<EvalType1, eval_type> || std::is_same_v<EvalType1, eval_deriv_type>);
1178 static_assert(
1179 std::is_same_v<EvalType2, eval_type> || std::is_same_v<EvalType2, eval_deriv_type>);
1180 ddc::DiscreteElement<bsplines_type1> jmin1;
1181 ddc::DiscreteElement<bsplines_type2> jmin2;
1182
1183 std::array<double, bsplines_type1::degree() + 1> vals1_ptr;
1184 Kokkos::mdspan<double, Kokkos::extents<std::size_t, bsplines_type1::degree() + 1>> const
1185 vals1(vals1_ptr.data());
1186 std::array<double, bsplines_type2::degree() + 1> vals2_ptr;
1187 Kokkos::mdspan<double, Kokkos::extents<std::size_t, bsplines_type2::degree() + 1>> const
1188 vals2(vals2_ptr.data());
1189 ddc::Coordinate<continuous_dimension_type1> const coord_eval_interest1
1190 = ddc::select<continuous_dimension_type1>(coord_eval);
1191 ddc::Coordinate<continuous_dimension_type2> const coord_eval_interest2
1192 = ddc::select<continuous_dimension_type2>(coord_eval);
1193
1194 if constexpr (std::is_same_v<EvalType1, eval_type>) {
1195 jmin1 = ddc::discrete_space<bsplines_type1>().eval_basis(vals1, coord_eval_interest1);
1196 } else if constexpr (std::is_same_v<EvalType1, eval_deriv_type>) {
1197 jmin1 = ddc::discrete_space<bsplines_type1>().eval_deriv(vals1, coord_eval_interest1);
1198 }
1199 if constexpr (std::is_same_v<EvalType2, eval_type>) {
1200 jmin2 = ddc::discrete_space<bsplines_type2>().eval_basis(vals2, coord_eval_interest2);
1201 } else if constexpr (std::is_same_v<EvalType2, eval_deriv_type>) {
1202 jmin2 = ddc::discrete_space<bsplines_type2>().eval_deriv(vals2, coord_eval_interest2);
1203 }
1204
1205 double y = 0.0;
1206 for (std::size_t i = 0; i < bsplines_type1::degree() + 1; ++i) {
1207 for (std::size_t j = 0; j < bsplines_type2::degree() + 1; ++j) {
1208 y += spline_coef(ddc::DiscreteElement<
1209 bsplines_type1,
1210 bsplines_type2>(jmin1 + i, jmin2 + j))
1211 * vals1[i] * vals2[j];
1212 }
1213 }
1214 return y;
1215 }
1216};
1217
1218} // namespace ddc
friend class DiscreteDomain
A class to evaluate, differentiate or integrate a 2D spline function.
void deriv_dim_1(ddc::ChunkSpan< double, batched_evaluation_domain_type, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< ddc::Coordinate< CoordsDims... > const, batched_evaluation_domain_type, Layout2, memory_space > const coords_eval, ddc::ChunkSpan< double const, batched_spline_domain_type, Layout3, memory_space > const spline_coef) const
Differentiate 2D spline function (described by its spline coefficients) on a mesh along first dimensi...
void deriv(ddc::ChunkSpan< double, batched_evaluation_domain_type, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< ddc::Coordinate< CoordsDims... > const, batched_evaluation_domain_type, Layout2, memory_space > const coords_eval, ddc::ChunkSpan< double const, batched_spline_domain_type, Layout3, memory_space > const spline_coef) const
Differentiate spline function (described by its spline coefficients) on a mesh along a specified dime...
void deriv_1_and_2(ddc::ChunkSpan< double, batched_evaluation_domain_type, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< double const, batched_spline_domain_type, Layout2, memory_space > const spline_coef) const
Cross-differentiate 2D spline function (described by its spline coefficients) on a mesh along dimensi...
KOKKOS_FUNCTION double deriv_dim_2(ddc::Coordinate< CoordsDims... > const &coord_eval, ddc::ChunkSpan< double const, spline_domain_type, Layout, memory_space > const spline_coef) const
Differentiate 2D spline function (described by its spline coefficients) at a given coordinate along s...
void deriv_dim_2(ddc::ChunkSpan< double, batched_evaluation_domain_type, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< ddc::Coordinate< CoordsDims... > const, batched_evaluation_domain_type, Layout2, memory_space > const coords_eval, ddc::ChunkSpan< double const, batched_spline_domain_type, Layout3, memory_space > const spline_coef) const
Differentiate 2D spline function (described by its spline coefficients) on a mesh along second dimens...
KOKKOS_FUNCTION double deriv_1_and_2(ddc::Coordinate< CoordsDims... > const &coord_eval, ddc::ChunkSpan< double const, spline_domain_type, Layout, memory_space > const spline_coef) const
Cross-differentiate 2D spline function (described by its spline coefficients) at a given coordinate.
SplineEvaluator2D(SplineEvaluator2D const &x)=default
Copy-constructs.
SplineEvaluator2D & operator=(SplineEvaluator2D const &x)=default
Copy-assigns.
void deriv2(ddc::ChunkSpan< double, batched_evaluation_domain_type, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< ddc::Coordinate< CoordsDims... > const, batched_evaluation_domain_type, Layout2, memory_space > const coords_eval, ddc::ChunkSpan< double const, batched_spline_domain_type, Layout3, memory_space > const spline_coef) const
Double-differentiate 2D spline function (described by its spline coefficients) on a mesh along specif...
SplineEvaluator2D(LowerExtrapolationRule1 const &lower_extrap_rule1, UpperExtrapolationRule1 const &upper_extrap_rule1, LowerExtrapolationRule2 const &lower_extrap_rule2, UpperExtrapolationRule2 const &upper_extrap_rule2)
Build a SplineEvaluator2D acting on batched_spline_domain.
void deriv_1_and_2(ddc::ChunkSpan< double, batched_evaluation_domain_type, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< ddc::Coordinate< CoordsDims... > const, batched_evaluation_domain_type, Layout2, memory_space > const coords_eval, ddc::ChunkSpan< double const, batched_spline_domain_type, Layout3, memory_space > const spline_coef) const
Cross-differentiate 2D spline function (described by its spline coefficients) on a mesh along dimensi...
void deriv_dim_2(ddc::ChunkSpan< double, batched_evaluation_domain_type, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< double const, batched_spline_domain_type, Layout2, memory_space > const spline_coef) const
Differentiate 2D spline function (described by its spline coefficients) on a mesh along second dimens...
KOKKOS_FUNCTION double deriv(ddc::Coordinate< CoordsDims... > const &coord_eval, ddc::ChunkSpan< double const, spline_domain_type, Layout, memory_space > const spline_coef) const
Differentiate 2D spline function (described by its spline coefficients) at a given coordinate along a...
~SplineEvaluator2D()=default
Destructs.
void deriv_dim_1(ddc::ChunkSpan< double, batched_evaluation_domain_type, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< double const, batched_spline_domain_type, Layout2, memory_space > const spline_coef) const
Differentiate 2D spline function (described by its spline coefficients) on a mesh along first dimensi...
void deriv2(ddc::ChunkSpan< double, batched_evaluation_domain_type, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< double const, batched_spline_domain_type, Layout2, memory_space > const spline_coef) const
Double-differentiate 2D spline function (described by its spline coefficients) on a mesh along specif...
void operator()(ddc::ChunkSpan< double, batched_evaluation_domain_type, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< double const, batched_spline_domain_type, Layout2, memory_space > const spline_coef) const
Evaluate 2D spline function (described by its spline coefficients) on a mesh.
upper_extrapolation_rule_2_type upper_extrapolation_rule_dim_2() const
Get the upper extrapolation rule along the second dimension.
lower_extrapolation_rule_2_type lower_extrapolation_rule_dim_2() const
Get the lower extrapolation rule along the second dimension.
void deriv(ddc::ChunkSpan< double, batched_evaluation_domain_type, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< double const, batched_spline_domain_type, Layout2, memory_space > const spline_coef) const
Differentiate spline function (described by its spline coefficients) on a mesh along a specified dime...
upper_extrapolation_rule_1_type upper_extrapolation_rule_dim_1() const
Get the upper extrapolation rule along the first dimension.
void integrate(ddc::ChunkSpan< double, batch_domain_type, Layout1, memory_space > const integrals, ddc::ChunkSpan< double const, batched_spline_domain_type, Layout2, memory_space > const spline_coef) const
Perform batched 2D integrations of a spline function (described by its spline coefficients) along the...
void operator()(ddc::ChunkSpan< double, batched_evaluation_domain_type, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< ddc::Coordinate< CoordsDims... > const, batched_evaluation_domain_type, Layout2, memory_space > const coords_eval, ddc::ChunkSpan< double const, batched_spline_domain_type, Layout3, memory_space > const spline_coef) const
Evaluate 2D spline function (described by its spline coefficients) on a mesh.
SplineEvaluator2D(SplineEvaluator2D &&x)=default
Move-constructs.
KOKKOS_FUNCTION double deriv2(ddc::Coordinate< CoordsDims... > const &coord_eval, ddc::ChunkSpan< double const, spline_domain_type, Layout, memory_space > const spline_coef) const
Double-differentiate 2D spline function (described by its spline coefficients) at a given coordinate ...
lower_extrapolation_rule_1_type lower_extrapolation_rule_dim_1() const
Get the lower extrapolation rule along the first dimension.
KOKKOS_FUNCTION double operator()(ddc::Coordinate< CoordsDims... > const &coord_eval, ddc::ChunkSpan< double const, spline_domain_type, Layout, memory_space > const spline_coef) const
Evaluate 2D spline function (described by its spline coefficients) at a given coordinate.
KOKKOS_FUNCTION double deriv_dim_1(ddc::Coordinate< CoordsDims... > const &coord_eval, ddc::ChunkSpan< double const, spline_domain_type, Layout, memory_space > const spline_coef) const
Differentiate 2D spline function (described by its spline coefficients) at a given coordinate along f...
SplineEvaluator2D & operator=(SplineEvaluator2D &&x)=default
Move-assigns.
The top-level namespace of DDC.