DDC 0.10.0
Loading...
Searching...
No Matches
spline_traits.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 <type_traits>
8
14
15namespace ddc {
16
17template <class T>
18struct is_spline_builder : std::false_type
19{
20};
21
22template <
23 class ExecSpace,
24 class MemorySpace,
25 class BSplines,
26 class InterpolationDDim,
27 ddc::BoundCond BcLower,
28 ddc::BoundCond BcUpper,
29 SplineSolver Solver>
31 ExecSpace,
32 MemorySpace,
33 BSplines,
34 InterpolationDDim,
35 BcLower,
36 BcUpper,
37 Solver>> : std::true_type
38{
39};
40
41/**
42 * @brief A helper to check if T is a SplineBuilder
43 * @tparam T The type to be checked if is a SplineBuilder
44 */
45template <class T>
46inline constexpr bool is_spline_builder_v = is_spline_builder<T>::value;
47
48template <class T>
49struct is_spline_builder2d : std::false_type
50{
51};
52
53template <
54 class ExecSpace,
55 class MemorySpace,
56 class BSpline1,
57 class BSpline2,
58 class DDimI1,
59 class DDimI2,
60 ddc::BoundCond BcLower1,
61 ddc::BoundCond BcUpper1,
62 ddc::BoundCond BcLower2,
63 ddc::BoundCond BcUpper2,
64 ddc::SplineSolver Solver>
66 ExecSpace,
67 MemorySpace,
68 BSpline1,
69 BSpline2,
70 DDimI1,
71 DDimI2,
72 BcLower1,
73 BcUpper1,
74 BcLower2,
75 BcUpper2,
76 Solver>> : std::true_type
77{
78};
79
80/**
81 * @brief A helper to check if T is a SplineBuilder2D
82 * @tparam T The type to be checked if is a SplineBuilder2D
83 */
84template <class T>
85inline constexpr bool is_spline_builder2d_v = is_spline_builder2d<T>::value;
86
87template <class T>
88struct is_spline_evaluator : std::false_type
89{
90};
91
92template <
93 class ExecSpace,
94 class MemorySpace,
95 class BSplines,
96 class EvaluationDDim,
97 class LowerExtrapolationRule,
98 class UpperExtrapolationRule>
100 ExecSpace,
101 MemorySpace,
102 BSplines,
103 EvaluationDDim,
104 LowerExtrapolationRule,
105 UpperExtrapolationRule>> : std::true_type
106{
107};
108
109/**
110 * @brief A helper to check if T is a SplineEvaluator
111 * @tparam T The type to be checked if is a SplineEvaluator
112 */
113template <class T>
114inline constexpr bool is_spline_evaluator_v = is_spline_evaluator<T>::value;
115
116template <class T>
117struct is_spline_evaluator2d : std::false_type
118{
119};
120
121template <
122 class ExecSpace,
123 class MemorySpace,
124 class BSpline1,
125 class BSpline2,
126 class EvaluationDDim1,
127 class EvaluationDDim2,
128 class LowerExtrapolationRule1,
129 class UpperExtrapolationRule1,
130 class LowerExtrapolationRule2,
131 class UpperExtrapolationRule2>
133 ExecSpace,
134 MemorySpace,
135 BSpline1,
136 BSpline2,
137 EvaluationDDim1,
138 EvaluationDDim2,
139 LowerExtrapolationRule1,
140 UpperExtrapolationRule1,
141 LowerExtrapolationRule2,
142 UpperExtrapolationRule2>> : std::true_type
143{
144};
145
146/**
147 * @brief A helper to check if T is a SplineEvaluator2D
148 * @tparam T The type to be checked if is a SplineEvaluator2D
149 */
150template <class T>
151inline constexpr bool is_spline_evaluator2d_v = is_spline_evaluator2d<T>::value;
152
153template <class Builder, class Evaluator>
154struct is_evaluator_admissible : std::false_type
155{
156};
157
158template <
159 class ExecSpace,
160 class MemorySpace,
161 class BSplines,
162 class InterpolationDDim,
163 ddc::BoundCond BcLower,
164 ddc::BoundCond BcUpper,
165 SplineSolver Solver,
166 class LowerExtrapolationRule,
167 class UpperExtrapolationRule>
170 ExecSpace,
171 MemorySpace,
172 BSplines,
173 InterpolationDDim,
174 BcLower,
175 BcUpper,
176 Solver>,
178 ExecSpace,
179 MemorySpace,
180 BSplines,
181 InterpolationDDim,
182 LowerExtrapolationRule,
183 UpperExtrapolationRule>> : std::true_type
184{
185};
186
187template <
188 class ExecSpace,
189 class MemorySpace,
190 class BSplines1,
191 class BSplines2,
192 class DDimI1,
193 class DDimI2,
194 ddc::BoundCond BcLower1,
195 ddc::BoundCond BcUpper1,
196 ddc::BoundCond BcLower2,
197 ddc::BoundCond BcUpper2,
198 SplineSolver Solver,
199 class LowerExtrapolationRule1,
200 class UpperExtrapolationRule1,
201 class LowerExtrapolationRule2,
202 class UpperExtrapolationRule2>
205 ExecSpace,
206 MemorySpace,
207 BSplines1,
208 BSplines2,
209 DDimI1,
210 DDimI2,
211 BcLower1,
212 BcUpper1,
213 BcLower2,
214 BcUpper2,
215 Solver>,
217 ExecSpace,
218 MemorySpace,
219 BSplines1,
220 BSplines2,
221 DDimI1,
222 DDimI2,
223 LowerExtrapolationRule1,
224 UpperExtrapolationRule1,
225 LowerExtrapolationRule2,
226 UpperExtrapolationRule2>> : std::true_type
227{
228};
229
230/**
231 * @brief A helper to check if SplineEvaluator is admissible for SplineBuilder
232 * @tparam Builder The builder type to be checked if it is admissible for Evaluator
233 * @tparam Evaluator The evaluator type to be checked if it is admissible for Builder
234 */
235template <class Builder, class Evaluator>
236inline constexpr bool is_evaluator_admissible_v
238
239} // namespace ddc
friend class ChunkSpan
friend class Chunk
Definition chunk.hpp:83
friend class DiscreteDomain
KOKKOS_DEFAULTED_FUNCTION constexpr DiscreteElement()=default
KOKKOS_FUNCTION constexpr bool operator!=(DiscreteVector< OTags... > const &rhs) const noexcept
A class which provides helper functions to initialise the Greville points from a B-Spline definition.
static ddc::DiscreteDomain< Sampling > get_domain()
Get the domain which gives us access to all of the Greville points.
static auto get_sampling()
Get the UniformPointSampling defining the Greville points.
Helper class for the initialisation of the mesh of interpolation points.
static auto get_sampling()
Get the sampling of interpolation points.
static ddc::DiscreteDomain< Sampling > get_domain()
Get the domain which can be used to access the interpolation points in the sampling.
Storage class of the static attributes of the discrete dimension.
Impl & operator=(Impl &&x)=default
Move-assigns.
Impl(RandomIt breaks_begin, RandomIt breaks_end)
Constructs an Impl by iterating over a range of break points from begin to end.
KOKKOS_INLINE_FUNCTION ddc::Coordinate< CDim > rmin() const noexcept
Returns the coordinate of the first break point of the domain on which the B-splines are defined.
Impl(std::vector< ddc::Coordinate< CDim > > const &breaks)
Constructs an Impl using a std::vector.
KOKKOS_INLINE_FUNCTION discrete_element_type eval_basis(DSpan1D values, ddc::Coordinate< CDim > const &x) const
Evaluates non-zero B-splines at a given coordinate.
KOKKOS_INLINE_FUNCTION std::size_t size() const noexcept
Returns the number of elements necessary to construct a spline representation of a function.
Impl(Impl< DDim, OriginMemorySpace > const &impl)
Copy-constructs from another Impl with a different Kokkos memory space.
~Impl()=default
Destructs.
KOKKOS_INLINE_FUNCTION ddc::DiscreteDomain< knot_discrete_dimension_type > break_point_domain() const
Returns the discrete domain which describes the break points.
KOKKOS_INLINE_FUNCTION ddc::DiscreteElement< knot_discrete_dimension_type > get_last_support_knot(discrete_element_type const &ix) const
Returns the coordinate of the last support knot associated to a DiscreteElement identifying a B-splin...
Impl(Impl &&x)=default
Move-constructs.
Impl(std::initializer_list< ddc::Coordinate< CDim > > breaks)
Constructs an Impl using a brace-list, i.e.
KOKKOS_INLINE_FUNCTION discrete_element_type eval_basis_and_n_derivs(ddc::DSpan2D derivs, ddc::Coordinate< CDim > const &x, std::size_t n) const
Evaluates non-zero B-spline values and derivatives at a given coordinate.
KOKKOS_INLINE_FUNCTION std::size_t ncells() const noexcept
Returns the number of cells over which the B-splines are defined.
KOKKOS_INLINE_FUNCTION discrete_domain_type full_domain() const
Returns the discrete domain including eventual additional B-splines in the periodic case.
KOKKOS_INLINE_FUNCTION ddc::DiscreteElement< knot_discrete_dimension_type > get_first_support_knot(discrete_element_type const &ix) const
Returns the coordinate of the first support knot associated to a DiscreteElement identifying a B-spli...
KOKKOS_INLINE_FUNCTION std::size_t npoints() const noexcept
The number of break points.
KOKKOS_INLINE_FUNCTION std::size_t nbasis() const noexcept
Returns the number of basis functions.
Impl(Impl const &x)=default
Copy-constructs.
KOKKOS_INLINE_FUNCTION discrete_element_type eval_deriv(DSpan1D derivs, ddc::Coordinate< CDim > const &x) const
Evaluates non-zero B-spline derivatives at a given coordinate.
KOKKOS_INLINE_FUNCTION double length() const noexcept
Returns the length of the domain.
Impl & operator=(Impl const &x)=default
Copy-assigns.
KOKKOS_INLINE_FUNCTION ddc::Coordinate< CDim > rmax() const noexcept
Returns the coordinate of the last break point of the domain on which the B-splines are defined.
The type of a non-uniform 1D spline basis (B-spline).
static constexpr std::size_t degree() noexcept
The degree of B-splines.
static constexpr bool is_periodic() noexcept
Indicates if the B-splines are periodic or not.
static constexpr bool is_uniform() noexcept
Indicates if the B-splines are uniform or not (this is not the case here).
NonUniformPointSampling models a non-uniform discretization of the CDim segment .
A class for creating a 2D spline approximation of a function.
void operator()(ddc::ChunkSpan< double, batched_spline_domain_type< BatchedInterpolationDDom >, Layout, memory_space > spline, ddc::ChunkSpan< double const, BatchedInterpolationDDom, Layout, memory_space > vals, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type1< BatchedInterpolationDDom >, Layout, memory_space > > derivs_min1=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type1< BatchedInterpolationDDom >, Layout, memory_space > > derivs_max1=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type2< BatchedInterpolationDDom >, Layout, memory_space > > derivs_min2=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type2< BatchedInterpolationDDom >, Layout, memory_space > > derivs_max2=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_min1_min2=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_max1_min2=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_min1_max2=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_max1_max2=std::nullopt) const
Compute a 2D spline approximation of a function.
ddc::DiscreteDomain< bsplines_type1, bsplines_type2 > spline_domain() const noexcept
Get the 2D domain on which spline coefficients are defined.
BatchedInterpolationDDom batched_interpolation_domain(BatchedInterpolationDDom const &batched_interpolation_domain) const noexcept
Get the whole domain representing interpolation points.
SplineBuilder2D(SplineBuilder2D &&x)=default
Move-constructs.
SplineBuilder2D & operator=(SplineBuilder2D &&x)=default
Move-assigns.
SplineBuilder2D(interpolation_domain_type const &interpolation_domain, std::optional< std::size_t > cols_per_chunk=std::nullopt, std::optional< unsigned int > preconditioner_max_block_size=std::nullopt)
Build a SplineBuilder2D acting on interpolation_domain.
batched_spline_domain_type< BatchedInterpolationDDom > batched_spline_domain(BatchedInterpolationDDom const &batched_interpolation_domain) const noexcept
Get the whole domain on which spline coefficients are defined.
batch_domain_type< BatchedInterpolationDDom > batch_domain(BatchedInterpolationDDom const &batched_interpolation_domain) const noexcept
Get the batch domain.
~SplineBuilder2D()=default
Destructs.
SplineBuilder2D(SplineBuilder2D const &x)=delete
Copy-constructor is deleted.
interpolation_domain_type interpolation_domain() const noexcept
Get the domain for the 2D interpolation mesh used by this class.
SplineBuilder2D(BatchedInterpolationDDom const &batched_interpolation_domain, std::optional< std::size_t > cols_per_chunk=std::nullopt, std::optional< unsigned int > preconditioner_max_block_size=std::nullopt)
Build a SplineBuilder2D acting on the interpolation domain contained in batched_interpolation_domain.
SplineBuilder2D & operator=(SplineBuilder2D const &x)=delete
Copy-assignment is deleted.
A class for creating a 3D spline approximation of a function.
SplineBuilder3D(SplineBuilder3D const &x)=delete
Copy-constructor is deleted.
void operator()(ddc::ChunkSpan< double, batched_spline_domain_type< BatchedInterpolationDDom >, Layout, memory_space > spline, ddc::ChunkSpan< double const, BatchedInterpolationDDom, Layout, memory_space > vals, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type1< BatchedInterpolationDDom >, Layout, memory_space > > derivs_min1=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type1< BatchedInterpolationDDom >, Layout, memory_space > > derivs_max1=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type2< BatchedInterpolationDDom >, Layout, memory_space > > derivs_min2=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type2< BatchedInterpolationDDom >, Layout, memory_space > > derivs_max2=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type3< BatchedInterpolationDDom >, Layout, memory_space > > derivs_min3=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type3< BatchedInterpolationDDom >, Layout, memory_space > > derivs_max3=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type1_2< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_min1_min2=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type1_2< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_max1_min2=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type1_2< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_min1_max2=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type1_2< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_max1_max2=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type2_3< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_min2_min3=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type2_3< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_max2_min3=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type2_3< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_min2_max3=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type2_3< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_max2_max3=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type1_3< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_min1_min3=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type1_3< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_max1_min3=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type1_3< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_min1_max3=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type1_3< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_max1_max3=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_min1_min2_min3=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_max1_min2_min3=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_min1_max2_min3=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_max1_max2_min3=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_min1_min2_max3=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_max1_min2_max3=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_min1_max2_max3=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type< BatchedInterpolationDDom >, Layout, memory_space > > mixed_derivs_max1_max2_max3=std::nullopt) const
Compute a 3D spline approximation of a function.
~SplineBuilder3D()=default
Destructs.
batched_spline_domain_type< BatchedInterpolationDDom > batched_spline_domain(BatchedInterpolationDDom const &batched_interpolation_domain) const noexcept
Get the whole domain on which spline coefficients are defined.
SplineBuilder3D(BatchedInterpolationDDom const &batched_interpolation_domain, std::optional< std::size_t > cols_per_chunk=std::nullopt, std::optional< unsigned int > preconditioner_max_block_size=std::nullopt)
Build a SplineBuilder3D acting on the interpolation domain contained in batched_interpolation_domain.
ddc::DiscreteDomain< bsplines_type1, bsplines_type2, bsplines_type3 > spline_domain() const noexcept
Get the 3D domain on which spline coefficients are defined.
batch_domain_type< BatchedInterpolationDDom > batch_domain(BatchedInterpolationDDom const &batched_interpolation_domain) const noexcept
Get the batch domain.
interpolation_domain_type interpolation_domain() const noexcept
Get the domain for the 3D interpolation mesh used by this class.
SplineBuilder3D & operator=(SplineBuilder3D const &x)=delete
Copy-assignment is deleted.
SplineBuilder3D(interpolation_domain_type const &interpolation_domain, std::optional< std::size_t > cols_per_chunk=std::nullopt, std::optional< unsigned int > preconditioner_max_block_size=std::nullopt)
Build a SplineBuilder3D acting on interpolation_domain.
SplineBuilder3D(SplineBuilder3D &&x)=default
Move-constructs.
BatchedInterpolationDDom batched_interpolation_domain(BatchedInterpolationDDom const &batched_interpolation_domain) const noexcept
Get the whole domain representing interpolation points.
SplineBuilder3D & operator=(SplineBuilder3D &&x)=default
Move-assigns.
A class for creating a spline approximation of a function.
batched_derivs_domain_type< BatchedInterpolationDDom > batched_derivs_xmax_domain(BatchedInterpolationDDom const &batched_interpolation_domain) const noexcept
Get the whole domain on which derivatives on upper boundary are defined.
static constexpr SplineSolver s_spline_solver
The SplineSolver giving the backend used to perform the spline approximation.
batch_domain_type< BatchedInterpolationDDom > batch_domain(BatchedInterpolationDDom const &batched_interpolation_domain) const noexcept
Get the batch domain.
std::tuple< ddc::Chunk< double, ddc::DiscreteDomain< ddc::Deriv< typename InterpolationDDim::continuous_dimension_type > >, ddc::KokkosAllocator< double, OutMemorySpace > >, ddc::Chunk< double, ddc::DiscreteDomain< InterpolationDDim >, ddc::KokkosAllocator< double, OutMemorySpace > >, ddc::Chunk< double, ddc::DiscreteDomain< ddc::Deriv< typename InterpolationDDim::continuous_dimension_type > >, ddc::KokkosAllocator< double, OutMemorySpace > > > quadrature_coefficients() const
Compute the quadrature coefficients associated to the b-splines used by this SplineBuilder.
SplineBuilder & operator=(SplineBuilder const &x)=delete
Copy-assignment is deleted.
SplineBuilder(SplineBuilder &&x)=default
Move-constructs.
SplineBuilder(interpolation_domain_type const &interpolation_domain, std::optional< std::size_t > cols_per_chunk=std::nullopt, std::optional< unsigned int > preconditioner_max_block_size=std::nullopt)
Build a SplineBuilder acting on interpolation_domain.
static constexpr ddc::BoundCond s_bc_xmin
The boundary condition implemented at the lower bound.
static constexpr int s_nbc_xmin
The number of equations defining the boundary condition at the lower bound.
SplineBuilder & operator=(SplineBuilder &&x)=default
Move-assigns.
BatchedInterpolationDDom batched_interpolation_domain(BatchedInterpolationDDom const &batched_interpolation_domain) const noexcept
Get the whole domain representing interpolation points.
SplineBuilder(SplineBuilder const &x)=delete
Copy-constructor is deleted.
SplineBuilder(BatchedInterpolationDDom const &batched_interpolation_domain, std::optional< std::size_t > cols_per_chunk=std::nullopt, std::optional< unsigned int > preconditioner_max_block_size=std::nullopt)
Build a SplineBuilder acting on the interpolation domain contained by batched_interpolation_domain.
static constexpr bool s_odd
Indicates if the degree of the splines is odd or even.
interpolation_domain_type interpolation_domain() const noexcept
Get the domain for the 1D interpolation mesh used by this class.
batched_derivs_domain_type< BatchedInterpolationDDom > batched_derivs_xmin_domain(BatchedInterpolationDDom const &batched_interpolation_domain) const noexcept
Get the whole domain on which derivatives on lower boundary are defined.
static constexpr ddc::BoundCond s_bc_xmax
The boundary condition implemented at the upper bound.
static constexpr int s_nbc_xmax
The number of equations defining the boundary condition at the upper bound.
void operator()(ddc::ChunkSpan< double, batched_spline_domain_type< BatchedInterpolationDDom >, Layout, memory_space > spline, ddc::ChunkSpan< double const, BatchedInterpolationDDom, Layout, memory_space > vals, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type< BatchedInterpolationDDom >, Layout, memory_space > > derivs_xmin=std::nullopt, std::optional< ddc::ChunkSpan< double const, batched_derivs_domain_type< BatchedInterpolationDDom >, Layout, memory_space > > derivs_xmax=std::nullopt) const
Compute a spline approximation of a function.
batched_spline_domain_type< BatchedInterpolationDDom > batched_spline_domain(BatchedInterpolationDDom const &batched_interpolation_domain) const noexcept
Get the whole domain on which spline coefficients are defined.
ddc::DiscreteDomain< bsplines_type > spline_domain() const noexcept
Get the 1D domain on which spline coefficients are defined.
~SplineBuilder()=default
Destructs.
A class to evaluate, differentiate or integrate a 2D spline function.
void operator()(ddc::ChunkSpan< double, BatchedInterpolationDDom, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< double const, batched_spline_domain_type< BatchedInterpolationDDom >, Layout2, memory_space > const spline_coef) const
Evaluate 2D spline function (described by its spline coefficients) on a mesh.
SplineEvaluator2D(SplineEvaluator2D &&x)=default
Move-constructs.
lower_extrapolation_rule_1_type lower_extrapolation_rule_dim_1() const
Get the lower extrapolation rule along the first dimension.
SplineEvaluator2D & operator=(SplineEvaluator2D const &x)=default
Copy-assigns.
void deriv(DElem const &deriv_order, ddc::ChunkSpan< double, BatchedInterpolationDDom, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< ddc::Coordinate< CoordsDims... > const, BatchedInterpolationDDom, Layout2, memory_space > const coords_eval, ddc::ChunkSpan< double const, batched_spline_domain_type< BatchedInterpolationDDom >, Layout3, memory_space > const spline_coef) const
Differentiate 2D spline function (described by its spline coefficients) on a mesh along the dimension...
SplineEvaluator2D(SplineEvaluator2D const &x)=default
Copy-constructs.
~SplineEvaluator2D()=default
Destructs.
KOKKOS_FUNCTION double deriv(DElem const &deriv_order, 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 t...
void deriv(DElem const &deriv_order, ddc::ChunkSpan< double, BatchedInterpolationDDom, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< double const, batched_spline_domain_type< BatchedInterpolationDDom >, Layout2, memory_space > const spline_coef) const
Differentiate 2D spline function (described by its spline coefficients) on a mesh along the dimension...
void integrate(ddc::ChunkSpan< double, BatchedDDom, Layout1, memory_space > const integrals, ddc::ChunkSpan< double const, BatchedSplineDDom, Layout2, memory_space > const spline_coef) const
Perform batched 2D integrations of a spline function (described by its spline coefficients) along the...
upper_extrapolation_rule_2_type upper_extrapolation_rule_dim_2() const
Get the upper extrapolation rule along the second dimension.
upper_extrapolation_rule_1_type upper_extrapolation_rule_dim_1() const
Get the upper extrapolation rule along the first dimension.
lower_extrapolation_rule_2_type lower_extrapolation_rule_dim_2() const
Get the lower extrapolation rule along the second dimension.
SplineEvaluator2D & operator=(SplineEvaluator2D &&x)=default
Move-assigns.
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.
void operator()(ddc::ChunkSpan< double, BatchedInterpolationDDom, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< ddc::Coordinate< CoordsDims... > const, BatchedInterpolationDDom, Layout2, memory_space > const coords_eval, ddc::ChunkSpan< double const, batched_spline_domain_type< BatchedInterpolationDDom >, Layout3, memory_space > const spline_coef) const
Evaluate 2D spline function (described by its spline coefficients) on a mesh.
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.
A class to evaluate, differentiate or integrate a 3D spline function.
upper_extrapolation_rule_1_type upper_extrapolation_rule_dim_1() const
Get the upper extrapolation rule along the first dimension.
void deriv(DElem const &deriv_order, ddc::ChunkSpan< double, BatchedInterpolationDDom, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< double const, batched_spline_domain_type< BatchedInterpolationDDom >, Layout2, memory_space > const spline_coef) const
Differentiate 3D spline function (described by its spline coefficients) on a mesh along the dimension...
void operator()(ddc::ChunkSpan< double, BatchedInterpolationDDom, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< ddc::Coordinate< CoordsDims... > const, BatchedInterpolationDDom, Layout2, memory_space > const coords_eval, ddc::ChunkSpan< double const, batched_spline_domain_type< BatchedInterpolationDDom >, Layout3, memory_space > const spline_coef) const
Evaluate 3D spline function (described by its spline coefficients) on a mesh.
SplineEvaluator3D & operator=(SplineEvaluator3D &&x)=default
Move-assigns.
upper_extrapolation_rule_3_type upper_extrapolation_rule_dim_3() const
Get the upper extrapolation rule along the third dimension.
SplineEvaluator3D(LowerExtrapolationRule1 const &lower_extrap_rule1, UpperExtrapolationRule1 const &upper_extrap_rule1, LowerExtrapolationRule2 const &lower_extrap_rule2, UpperExtrapolationRule2 const &upper_extrap_rule2, LowerExtrapolationRule3 const &lower_extrap_rule3, UpperExtrapolationRule3 const &upper_extrap_rule3)
Build a SplineEvaluator3D acting on batched_spline_domain.
~SplineEvaluator3D()=default
Destructs.
lower_extrapolation_rule_1_type lower_extrapolation_rule_dim_1() const
Get the lower extrapolation rule along the first dimension.
upper_extrapolation_rule_2_type upper_extrapolation_rule_dim_2() const
Get the upper extrapolation rule along the second dimension.
void deriv(DElem const &deriv_order, ddc::ChunkSpan< double, BatchedInterpolationDDom, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< ddc::Coordinate< CoordsDims... > const, BatchedInterpolationDDom, Layout2, memory_space > const coords_eval, ddc::ChunkSpan< double const, batched_spline_domain_type< BatchedInterpolationDDom >, Layout3, memory_space > const spline_coef) const
Differentiate 3D spline function (described by its spline coefficients) on a mesh along the dimension...
KOKKOS_FUNCTION double deriv(DElem const &deriv_order, ddc::Coordinate< CoordsDims... > const &coord_eval, ddc::ChunkSpan< double const, spline_domain_type, Layout, memory_space > const spline_coef) const
Differentiate 3D spline function (described by its spline coefficients) at a given coordinate along t...
SplineEvaluator3D(SplineEvaluator3D const &x)=default
Copy-constructs.
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 3D spline function (described by its spline coefficients) at a given coordinate.
lower_extrapolation_rule_3_type lower_extrapolation_rule_dim_3() const
Get the lower extrapolation rule along the third dimension.
SplineEvaluator3D(SplineEvaluator3D &&x)=default
Move-constructs.
void integrate(ddc::ChunkSpan< double, BatchedDDom, Layout1, memory_space > const integrals, ddc::ChunkSpan< double const, BatchedSplineDDom, Layout2, memory_space > const spline_coef) const
Perform batched 3D integrations of a spline function (described by its spline coefficients) along the...
lower_extrapolation_rule_2_type lower_extrapolation_rule_dim_2() const
Get the lower extrapolation rule along the second dimension.
SplineEvaluator3D & operator=(SplineEvaluator3D const &x)=default
Copy-assigns.
void operator()(ddc::ChunkSpan< double, BatchedInterpolationDDom, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< double const, batched_spline_domain_type< BatchedInterpolationDDom >, Layout2, memory_space > const spline_coef) const
Evaluate 3D spline function (described by its spline coefficients) on a mesh.
A class to evaluate, differentiate or integrate a spline function.
void operator()(ddc::ChunkSpan< double, BatchedInterpolationDDom, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< ddc::Coordinate< CoordsDims... > const, BatchedInterpolationDDom, Layout2, memory_space > const coords_eval, ddc::ChunkSpan< double const, batched_spline_domain_type< BatchedInterpolationDDom >, Layout3, memory_space > const spline_coef) const
Evaluate spline function (described by its spline coefficients) on a mesh.
upper_extrapolation_rule_type upper_extrapolation_rule() const
Get the upper extrapolation rule.
SplineEvaluator & operator=(SplineEvaluator const &x)=default
Copy-assigns.
SplineEvaluator & operator=(SplineEvaluator &&x)=default
Move-assigns.
SplineEvaluator(LowerExtrapolationRule const &lower_extrap_rule, UpperExtrapolationRule const &upper_extrap_rule)
Build a SplineEvaluator acting on batched_spline_domain.
lower_extrapolation_rule_type lower_extrapolation_rule() const
Get the lower extrapolation rule.
void deriv(DElem const &deriv_order, ddc::ChunkSpan< double, BatchedInterpolationDDom, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< double const, batched_spline_domain_type< BatchedInterpolationDDom >, Layout2, memory_space > const spline_coef) const
Differentiate 1D spline function (described by its spline coefficients) on a mesh.
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 1D spline function (described by its spline coefficients) at a given coordinate.
KOKKOS_FUNCTION double deriv(DElem const &deriv_order, ddc::Coordinate< CoordsDims... > const &coord_eval, ddc::ChunkSpan< double const, spline_domain_type, Layout, memory_space > const spline_coef) const
Differentiate 1D spline function (described by its spline coefficients) at a given coordinate.
SplineEvaluator(SplineEvaluator const &x)=default
Copy-constructs.
SplineEvaluator(SplineEvaluator &&x)=default
Move-constructs.
void deriv(DElem const &deriv_order, ddc::ChunkSpan< double, BatchedInterpolationDDom, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< ddc::Coordinate< CoordsDims... > const, BatchedInterpolationDDom, Layout2, memory_space > const coords_eval, ddc::ChunkSpan< double const, batched_spline_domain_type< BatchedInterpolationDDom >, Layout3, memory_space > const spline_coef) const
Differentiate 1D spline function (described by its spline coefficients) on a mesh.
void integrate(ddc::ChunkSpan< double, BatchedDDom, Layout1, memory_space > const integrals, ddc::ChunkSpan< double const, BatchedSplineDDom, Layout2, memory_space > const spline_coef) const
Perform batched 1D integrations of a spline function (described by its spline coefficients) along the...
void operator()(ddc::ChunkSpan< double, BatchedInterpolationDDom, Layout1, memory_space > const spline_eval, ddc::ChunkSpan< double const, batched_spline_domain_type< BatchedInterpolationDDom >, Layout2, memory_space > const spline_coef) const
Evaluate a spline function (described by its spline coefficients) on a mesh.
~SplineEvaluator()=default
Destructs.
Storage class of the static attributes of the discrete dimension.
KOKKOS_INLINE_FUNCTION ddc::DiscreteElement< knot_discrete_dimension_type > get_last_support_knot(discrete_element_type const &ix) const
Returns the coordinate of the last support knot associated to a DiscreteElement identifying a B-splin...
Impl(ddc::Coordinate< CDim > rmin, ddc::Coordinate< CDim > rmax, std::size_t ncells)
Constructs a spline basis (B-splines) with n equidistant knots over .
KOKKOS_INLINE_FUNCTION ddc::Coordinate< CDim > rmax() const noexcept
Returns the coordinate of the upper bound of the domain on which the B-splines are defined.
KOKKOS_INLINE_FUNCTION discrete_element_type eval_basis(DSpan1D values, ddc::Coordinate< CDim > const &x) const
Evaluates non-zero B-splines at a given coordinate.
KOKKOS_INLINE_FUNCTION ddc::DiscreteDomain< knot_discrete_dimension_type > break_point_domain() const
Returns the discrete domain which describes the break points.
KOKKOS_INLINE_FUNCTION ddc::Coordinate< CDim > rmin() const noexcept
Returns the coordinate of the lower bound of the domain on which the B-splines are defined.
~Impl()=default
Destructs.
KOKKOS_INLINE_FUNCTION std::size_t nbasis() const noexcept
Returns the number of basis functions.
Impl(Impl const &x)=default
Copy-constructs.
KOKKOS_INLINE_FUNCTION std::size_t size() const noexcept
Returns the number of elements necessary to construct a spline representation of a function.
Impl(Impl &&x)=default
Move-constructs.
KOKKOS_INLINE_FUNCTION discrete_element_type eval_basis_and_n_derivs(ddc::DSpan2D derivs, ddc::Coordinate< CDim > const &x, std::size_t n) const
Evaluates non-zero B-spline values and derivatives at a given coordinate.
KOKKOS_INLINE_FUNCTION ddc::DiscreteElement< knot_discrete_dimension_type > get_first_support_knot(discrete_element_type const &ix) const
Returns the coordinate of the first support knot associated to a DiscreteElement identifying a B-spli...
KOKKOS_INLINE_FUNCTION double length() const noexcept
Returns the length of the domain.
Impl(Impl< DDim, OriginMemorySpace > const &impl)
Copy-constructs from another Impl with a different Kokkos memory space.
KOKKOS_INLINE_FUNCTION std::size_t ncells() const noexcept
Returns the number of cells over which the B-splines are defined.
KOKKOS_INLINE_FUNCTION discrete_element_type eval_deriv(DSpan1D derivs, ddc::Coordinate< CDim > const &x) const
Evaluates non-zero B-spline derivatives at a given coordinate.
Impl & operator=(Impl &&x)=default
Move-assigns.
KOKKOS_INLINE_FUNCTION discrete_domain_type full_domain() const
Returns the discrete domain including eventual additional B-splines in the periodic case.
Impl & operator=(Impl const &x)=default
Copy-assigns.
The type of a uniform 1D spline basis (B-spline).
static constexpr bool is_uniform() noexcept
Indicates if the B-splines are uniform or not (this is the case here).
static constexpr bool is_periodic() noexcept
Indicates if the B-splines are periodic or not.
static constexpr std::size_t degree() noexcept
The degree of B-splines.
UniformPointSampling models a uniform discretization of the provided continuous dimension.
#define DDC_BUILD_DEPRECATED_CODE
Definition config.hpp:7
The top-level namespace of DDC.
constexpr int n_boundary_equations(ddc::BoundCond const bc, std::size_t const degree)
Return the number of equations needed to describe a given boundary condition.
constexpr bool is_spline_builder_v
A helper to check if T is a SplineBuilder.
constexpr bool is_uniform_bsplines_v
Indicates if a tag corresponds to uniform B-splines or not.
BoundCond
An enum representing a spline boundary condition.
@ GREVILLE
Use Greville points instead of conditions on derivative for B-Spline interpolation.
@ HERMITE
Hermite boundary condition.
@ PERIODIC
Periodic boundary condition u(1)=u(n)
constexpr bool is_spline_evaluator2d_v
A helper to check if T is a SplineEvaluator2D.
ddc::ChunkSpan< double, ddc::DiscreteDomain< DDim >, Layout, MemorySpace > integrals(ExecSpace const &execution_space, ddc::ChunkSpan< double, ddc::DiscreteDomain< DDim >, Layout, MemorySpace > int_vals)
Compute the integrals of the B-splines.
SplineSolver
An enum determining the backend solver of a SplineBuilder or SplineBuilder2d.
@ LAPACK
Enum member to identify the LAPACK-based solver (direct method)
@ GINKGO
Enum member to identify the Ginkgo-based solver (iterative method)
constexpr bool is_spline_builder2d_v
A helper to check if T is a SplineBuilder2D.
constexpr bool is_spline_evaluator_v
A helper to check if T is a SplineEvaluator.
constexpr bool is_non_uniform_bsplines_v
Indicates if a tag corresponds to non-uniform B-splines or not.
constexpr bool is_evaluator_admissible_v
A helper to check if SplineEvaluator is admissible for SplineBuilder.
A templated struct representing a discrete dimension storing the derivatives of a function along a co...
Definition deriv.hpp:15
If the type DDim is a B-spline, defines type to the discrete dimension of the associated knots.
ConstantExtrapolationRule(ddc::Coordinate< DimI > eval_pos, ddc::Coordinate< DimNI > eval_pos_not_interest_min, ddc::Coordinate< DimNI > eval_pos_not_interest_max)
Instantiate a ConstantExtrapolationRule.
KOKKOS_FUNCTION double operator()(CoordType coord_extrap, ddc::ChunkSpan< double const, ddc::DiscreteDomain< BSplines1, BSplines2 >, Layout, MemorySpace > const spline_coef) const
Get the value of the function on B-splines at a coordinate outside the domain.
ConstantExtrapolationRule(ddc::Coordinate< DimI > eval_pos)
Instantiate a ConstantExtrapolationRule.
KOKKOS_FUNCTION double operator()(CoordType pos, ddc::ChunkSpan< double const, ddc::DiscreteDomain< BSplines >, Layout, MemorySpace > const spline_coef) const
Get the value of the function on B-splines at a coordinate outside the domain.
ConstantExtrapolationRule(ddc::Coordinate< DimI > eval_pos)
Instantiate a ConstantExtrapolationRule.
A functor describing a null extrapolation boundary value for 1D spline evaluator.
KOKKOS_FUNCTION double operator()(CoordType, ChunkSpan) const
Evaluates the spline at a coordinate outside of the domain.
KOKKOS_FUNCTION double operator()(CoordType, ChunkSpan) const