Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Provide a new initialization interface #719

Draft
wants to merge 1 commit into
base: main
Choose a base branch
from
Draft

Provide a new initialization interface #719

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
12 changes: 5 additions & 7 deletions benchmarks/splines.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -111,26 +111,24 @@ void characteristics_advection_unitary(benchmark::State& state)
std::ref(maxUsedMem));

if constexpr (!IsNonUniform) {
ddc::init_discrete_space<BSplinesX<
ddc::create_uniform_bsplines<BSplinesX<
IsNonUniform,
s_degree_x>>(ddc::Coordinate<X>(0.), ddc::Coordinate<X>(1.), nx);
} else {
std::vector<ddc::Coordinate<X>> breaks(nx + 1);
for (std::size_t i(0); i < nx + 1; ++i) {
breaks[i] = ddc::Coordinate<X>(static_cast<double>(i) / nx);
}
ddc::init_discrete_space<BSplinesX<IsNonUniform, s_degree_x>>(breaks);
ddc::create_non_uniform_bsplines<BSplinesX<IsNonUniform, s_degree_x>>(breaks);
}
ddc::init_discrete_space<DDimX<IsNonUniform, s_degree_x>>(
ddc::init_discrete_space_from_impl<DDimX<IsNonUniform, s_degree_x>>(
ddc::GrevilleInterpolationPoints<
BSplinesX<IsNonUniform, s_degree_x>,
ddc::BoundCond::PERIODIC,
ddc::BoundCond::PERIODIC>::
template get_sampling<DDimX<IsNonUniform, s_degree_x>>());
ddc::DiscreteDomain<DDimY> const y_domain = ddc::init_discrete_space<DDimY>(DDimY::init<DDimY>(
ddc::Coordinate<Y>(-1.),
ddc::Coordinate<Y>(1.),
ddc::DiscreteVector<DDimY>(ny)));
ddc::DiscreteDomain<DDimY> const y_domain = ddc::create_uniform_point_sampling<
DDimY>(ddc::Coordinate<Y>(-1.), ddc::Coordinate<Y>(1.), ddc::DiscreteVector<DDimY>(ny));

auto const x_domain = ddc::GrevilleInterpolationPoints<
BSplinesX<IsNonUniform, s_degree_x>,
Expand Down
2 changes: 1 addition & 1 deletion docs/first_steps.md
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -75,7 +75,7 @@ Their type is not specified because we use C++
all of the same type: `DiscreteDomain<DDimX>` that represents a set of contiguous points in the
discretization of `X`.

\ref ddc::UniformPointSampling::init_ghosted "init_ghosted" takes as parameters the coordinate of the first and last discretized points, the
\ref ddc::create_uniform_point_sampling_ghosted "init_ghosted" takes as parameters the coordinate of the first and last discretized points, the
number of discretized points in the domain and the number of additional points on each side of the
domain.
The fours `DiscreteDomain`s returned are:
Expand Down
15 changes: 6 additions & 9 deletions examples/characteristics_advection.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -133,17 +133,17 @@ int main(int argc, char** argv)

//! [X-global-domain]
// Initialization of the global domain in X
ddc::init_discrete_space<
ddc::create_uniform_bsplines<
BSplinesX>(ddc::Coordinate<X>(x_start), ddc::Coordinate<X>(x_end), nb_x_points);
ddc::init_discrete_space<DDimX>(GrevillePoints::get_sampling<DDimX>());
ddc::init_discrete_space_from_impl<DDimX>(GrevillePoints::get_sampling<DDimX>());

auto const x_domain = GrevillePoints::get_domain<DDimX>();
//! [X-global-domain]
// Initialization of the global domain in Y
auto const y_domain = ddc::init_discrete_space<DDimY>(DDimY::init<DDimY>(
ddc::DiscreteDomain<DDimY> const y_domain = ddc::create_uniform_point_sampling<DDimY>(
ddc::Coordinate<Y>(y_start),
ddc::Coordinate<Y>(y_end),
ddc::DiscreteVector<DDimY>(nb_y_points)));
ddc::DiscreteVector<DDimY>(nb_y_points));

//! [time-domains]

Expand All @@ -153,11 +153,8 @@ int main(int argc, char** argv)
// Initialization of the global domain in time:
// - the number of discrete time-points is equal to the number of
// steps + 1
ddc::DiscreteDomain<DDimT> const time_domain
= ddc::init_discrete_space<DDimT>(DDimT::init<DDimT>(
ddc::Coordinate<T>(start_time),
ddc::Coordinate<T>(end_time),
nb_time_steps + 1));
ddc::DiscreteDomain<DDimT> const time_domain = ddc::create_uniform_point_sampling<
DDimT>(ddc::Coordinate<T>(start_time), ddc::Coordinate<T>(end_time), nb_time_steps + 1);
//! [time-domains]

//! [data allocation]
Expand Down
15 changes: 6 additions & 9 deletions examples/heat_equation.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -122,11 +122,11 @@ int main(int argc, char** argv)
// Initialization of the global domain in X with gwx ghost points on
// each side
auto const [x_domain, ghosted_x_domain, x_pre_ghost, x_post_ghost]
= ddc::init_discrete_space<DDimX>(DDimX::init_ghosted<DDimX>(
= ddc::create_uniform_point_sampling_ghosted<DDimX>(
ddc::Coordinate<X>(x_start),
ddc::Coordinate<X>(x_end),
ddc::DiscreteVector<DDimX>(nb_x_points),
gwx));
gwx);
//! [X-global-domain]

//! [X-domains]
Expand All @@ -147,11 +147,11 @@ int main(int argc, char** argv)
// Initialization of the global domain in Y with gwy ghost points on
// each side
auto const [y_domain, ghosted_y_domain, y_pre_ghost, y_post_ghost]
= ddc::init_discrete_space<DDimY>(DDimY::init_ghosted<DDimY>(
= ddc::create_uniform_point_sampling_ghosted<DDimY>(
ddc::Coordinate<Y>(y_start),
ddc::Coordinate<Y>(y_end),
ddc::DiscreteVector<DDimY>(nb_y_points),
gwy));
gwy);

// our zone at the start of the domain that will be mirrored to the
// ghost
Expand Down Expand Up @@ -190,11 +190,8 @@ int main(int argc, char** argv)
// - the number of discrete time-points is equal to the number of
// steps + 1
// `init` takes required information to initialize the attributes of the dimension.
ddc::DiscreteDomain<DDimT> const time_domain
= ddc::init_discrete_space<DDimT>(DDimT::init<DDimT>(
ddc::Coordinate<T>(start_time),
ddc::Coordinate<T>(end_time),
nb_time_steps + 1));
ddc::DiscreteDomain<DDimT> const time_domain = ddc::create_uniform_point_sampling<
DDimT>(ddc::Coordinate<T>(start_time), ddc::Coordinate<T>(end_time), nb_time_steps + 1);
//! [time-domains]

//! [data allocation]
Expand Down
19 changes: 8 additions & 11 deletions examples/heat_equation_spectral.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -99,18 +99,18 @@ int main(int argc, char** argv)
std::ptrdiff_t const t_output_period = 10;

// Initialization of the global domain in X including periodic point to have correct step
auto const x_domain_with_periodic_point = ddc::init_discrete_space<DDimX>(DDimX::init<DDimX>(
auto const x_domain_with_periodic_point = ddc::create_uniform_point_sampling<DDimX>(
ddc::Coordinate<X>(x_start),
ddc::Coordinate<X>(x_end),
ddc::DiscreteVector<DDimX>(nb_x_points)));
ddc::DiscreteVector<DDimX>(nb_x_points));
ddc::DiscreteDomain<DDimX> const x_domain
= x_domain_with_periodic_point.remove_last(ddc::DiscreteVector<DDimX>(1));

// Initialization of the global domain in Y including periodic point to have correct step
auto const y_domain_with_periodic_point = ddc::init_discrete_space<DDimY>(DDimY::init<DDimY>(
auto const y_domain_with_periodic_point = ddc::create_uniform_point_sampling<DDimY>(
ddc::Coordinate<Y>(y_start),
ddc::Coordinate<Y>(y_end),
ddc::DiscreteVector<DDimY>(nb_y_points)));
ddc::DiscreteVector<DDimY>(nb_y_points));
ddc::DiscreteDomain<DDimY> const y_domain
= y_domain_with_periodic_point.remove_last(ddc::DiscreteVector<DDimY>(1));

Expand All @@ -123,11 +123,8 @@ int main(int argc, char** argv)
// Initialization of the global domain in time:
// - the number of discrete time-points is equal to the number of
// steps + 1
ddc::DiscreteDomain<DDimT> const time_domain
= ddc::init_discrete_space<DDimT>(DDimT::init<DDimT>(
ddc::Coordinate<T>(start_time),
ddc::Coordinate<T>(end_time),
nb_time_steps + 1));
ddc::DiscreteDomain<DDimT> const time_domain = ddc::create_uniform_point_sampling<
DDimT>(ddc::Coordinate<T>(start_time), ddc::Coordinate<T>(end_time), nb_time_steps + 1);

ddc::DiscreteDomain<DDimX, DDimY> const xy_domain(x_domain, y_domain);

Expand Down Expand Up @@ -156,8 +153,8 @@ int main(int argc, char** argv)
// time of the iteration where the last output happened
ddc::DiscreteElement<DDimT> last_output = time_domain.front();

ddc::init_discrete_space<DDimFx>(ddc::init_fourier_space<DDimFx>(x_domain));
ddc::init_discrete_space<DDimFy>(ddc::init_fourier_space<DDimFy>(y_domain));
ddc::create_fourier_space<DDimFx>(x_domain);
ddc::create_fourier_space<DDimFy>(y_domain);

ddc::DiscreteDomain<DDimFx, DDimFy> const k_mesh
= ddc::fourier_mesh<DDimFx, DDimFy>(xy_domain, false);
Expand Down
15 changes: 6 additions & 9 deletions examples/non_uniform_heat_equation.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -173,8 +173,8 @@ int main(int argc, char** argv)

//! [build-domains]
auto const [x_domain, ghosted_x_domain, x_pre_ghost, x_post_ghost]
= ddc::init_discrete_space<DDimX>(
DDimX::init_ghosted<DDimX>(x_domain_vect, x_pre_ghost_vect, x_post_ghost_vect));
= ddc::create_non_uniform_point_sampling_ghosted<
DDimX>(x_domain_vect, x_pre_ghost_vect, x_post_ghost_vect);
//! [build-domains]

ddc::DiscreteDomain<DDimX> const
Expand All @@ -193,8 +193,8 @@ int main(int argc, char** argv)

//! [build-Y-domain]
auto const [y_domain, ghosted_y_domain, y_pre_ghost, y_post_ghost]
= ddc::init_discrete_space<DDimY>(
DDimY::init_ghosted<DDimY>(y_domain_vect, y_pre_ghost_vect, y_post_ghost_vect));
= ddc::create_non_uniform_point_sampling_ghosted<
DDimY>(y_domain_vect, y_pre_ghost_vect, y_post_ghost_vect);
//! [build-Y-domain]

ddc::DiscreteDomain<DDimY> const
Expand Down Expand Up @@ -227,11 +227,8 @@ int main(int argc, char** argv)
ddc::DiscreteVector<DDimT> const nb_time_steps(
std::ceil((end_time - start_time) / max_dt) + .2);

ddc::DiscreteDomain<DDimT> const time_domain
= ddc::init_discrete_space<DDimT>(DDimT::init<DDimT>(
ddc::Coordinate<T>(start_time),
ddc::Coordinate<T>(end_time),
nb_time_steps + 1));
ddc::DiscreteDomain<DDimT> const time_domain = ddc::create_uniform_point_sampling<
DDimT>(ddc::Coordinate<T>(start_time), ddc::Coordinate<T>(end_time), nb_time_steps + 1);
//! [time-domain]

//! [data allocation]
Expand Down
15 changes: 6 additions & 9 deletions examples/uniform_heat_equation.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -96,11 +96,11 @@ int main(int argc, char** argv)

//! [X-global-domain]
auto const [x_domain, ghosted_x_domain, x_pre_ghost, x_post_ghost]
= ddc::init_discrete_space<DDimX>(DDimX::init_ghosted<DDimX>(
= ddc::create_uniform_point_sampling_ghosted<DDimX>(
ddc::Coordinate<X>(x_start),
ddc::Coordinate<X>(x_end),
ddc::DiscreteVector<DDimX>(nb_x_points),
gwx));
gwx);
//! [X-global-domain]

//! [X-domains]
Expand All @@ -114,11 +114,11 @@ int main(int argc, char** argv)
ddc::DiscreteVector<DDimY> const gwy(1);

auto const [y_domain, ghosted_y_domain, y_pre_ghost, y_post_ghost]
= ddc::init_discrete_space<DDimY>(DDimY::init_ghosted<DDimY>(
= ddc::create_uniform_point_sampling_ghosted<DDimY>(
ddc::Coordinate<Y>(y_start),
ddc::Coordinate<Y>(y_end),
ddc::DiscreteVector<DDimY>(nb_y_points),
gwy));
gwy);

ddc::DiscreteDomain<DDimY> const
y_post_mirror(y_post_ghost.front() - y_domain.extents(), y_post_ghost.extents());
Expand All @@ -139,11 +139,8 @@ int main(int argc, char** argv)
//! [time-domain]
ddc::DiscreteVector<DDimT> const nb_time_steps(std::ceil((end_time - start_time) / dt) + .2);

ddc::DiscreteDomain<DDimT> const time_domain
= ddc::init_discrete_space<DDimT>(DDimT::init<DDimT>(
ddc::Coordinate<T>(start_time),
ddc::Coordinate<T>(end_time),
nb_time_steps + 1));
ddc::DiscreteDomain<DDimT> const time_domain = ddc::create_uniform_point_sampling<
DDimT>(ddc::Coordinate<T>(start_time), ddc::Coordinate<T>(end_time), nb_time_steps + 1);
//! [time-domain]

//! [data allocation]
Expand Down
10 changes: 10 additions & 0 deletions include/ddc/discrete_space.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -216,6 +216,16 @@ std::tuple<Arg0, Arg1, Args...> init_discrete_space(std::tuple<DDimImpl, Arg0, A
return detail::extract_after(std::move(a), std::index_sequence_for<Arg0, Arg1, Args...>());
}

/** Initialize (emplace) a global singleton discrete space
*
* @param ddim_impl the impl discrete space object
*/
template <class DDim, class DDimImpl>
void init_discrete_space_from_impl(DDimImpl&& ddim_impl)
{
init_discrete_space<DDim>(std::forward<DDimImpl>(ddim_impl));
}

/**
* @tparam DDim a discrete dimension
* @return the discrete space instance associated with `DDim`.
Expand Down
8 changes: 8 additions & 0 deletions include/ddc/experimental/single_discretization.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -92,6 +92,14 @@ class SingleDiscretization : SingleDiscretizationBase
};
};

template <class DDim>
void create_single_discretization(Coordinate<typename DDim::continuous_dimension_type> origin)
{
init_discrete_space_from_impl<DDim>(
typename SingleDiscretization<typename DDim::continuous_dimension_type>::
template Impl<DDim, Kokkos::HostSpace>(origin));
}

template <class DDim>
struct is_single_discretization : public std::is_base_of<SingleDiscretizationBase, DDim>::type
{
Expand Down
41 changes: 41 additions & 0 deletions include/ddc/kernels/fft.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -300,6 +300,47 @@ typename DDimFx::template Impl<DDimFx, Kokkos::HostSpace> init_fourier_space(
return std::move(impl);
}

/**
* @brief Initialize a Fourier discrete dimension.
*
* Initialize the (1D) discrete space representing the Fourier discrete dimension associated
* to the (1D) mesh passed as argument. It is a N-periodic PeriodicSampling with a periodic window of width 2*pi/dx.
*
* This value comes from the Nyquist-Shannon theorem: the period of the spectral domain is N*dk = 2*pi/dx.
* Adding to this the relations dx = (xmax-xmin)/(N-1), and dk = (kmax-kmin)/(N-1), we get kmax-kmin = 2*pi*(N-1)^2/N/(xmax-xmin),
* which is used in the implementation (xmax, xmin, kmin and kmax are the centers of lower and upper cells inside a single period of the meshes).
*
* @tparam DDimFx A PeriodicSampling representing the Fourier discrete dimension.
* @tparam DDimX The type of the original discrete dimension.
*
* @param x_mesh The DiscreteDomain representing the (1D) original mesh.
*
* @see PeriodicSampling
*/
template <typename DDimFx, typename DDimX>
void create_fourier_space(ddc::DiscreteDomain<DDimX> x_mesh)
{
static_assert(
is_uniform_point_sampling_v<DDimX>,
"DDimX dimension must derive from UniformPointSampling");
static_assert(
is_periodic_sampling_v<DDimFx>,
"DDimFx dimension must derive from PeriodicSampling");
using CDimFx = typename DDimFx::continuous_dimension_type;
using CDimX = typename DDimX::continuous_dimension_type;
static_assert(
std::is_same_v<CDimFx, ddc::Fourier<CDimX>>,
"DDimX and DDimFx dimensions must be defined over the same continuous dimension");

DiscreteVectorElement const nx = get<DDimX>(x_mesh.extents());
double const lx = ddc::rlength(x_mesh);
create_periodic_sampling<DDimFx>(
ddc::Coordinate<CDimFx>(0),
ddc::Coordinate<CDimFx>(2 * (nx - 1) * (nx - 1) / (nx * lx) * Kokkos::numbers::pi),
ddc::DiscreteVector<DDimFx>(nx),
ddc::DiscreteVector<DDimFx>(nx));
}

/**
* @brief Get the Fourier mesh.
*
Expand Down
42 changes: 41 additions & 1 deletion include/ddc/kernels/splines/bsplines_non_uniform.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -385,6 +385,46 @@ class NonUniformBSplines : detail::NonUniformBSplinesBase
};
};

/** @brief Constructs a non uniform spline basis by iterating over a range of break points from begin to end.
* The associated static attributes are available after this call.
*
* The provided break points describe the separation between the cells on which the polynomials
* comprising a spline are defined. They are used to build a set of knots. There are 2*degree more
* knots than break points. In the non-periodic case the knots are defined as follows:
* \f$ k_i = b_0 \forall 0 \leq i < d \f$
* \f$ k_{i+d} = b_i \forall 0 \leq i < n_b \f$
* \f$ k_{i+d+n_b} = b_{n_b-1} \forall 0 \leq i < d \f$
* where \f$d\f$ is the degree of the polynomials, and \f$n_b\f$ is the number of break points in the input pair of iterators. And in the periodic case:
* \f$ k_i = b_{n_b-1-d+i} \forall 0 \leq i < d \f$
* \f$ k_{i+d} = b_i \forall 0 \leq i \leq n_b \f$
* \f$ k_{i+d+n_b} = b_{i+1} \forall 0 \leq i < d \f$
*
* This constructor makes the knots accessible via a DiscreteSpace.
* @tparam DDim The name of the discrete dimension
* @param breaks_begin The iterator which points at the beginning of the break points.
* @param breaks_end The iterator which points at the end of the break points.
*/
template <class DDim, class RandomIt>
void create_non_uniform_bsplines(RandomIt const breaks_begin, RandomIt const breaks_end)
{
ddc::init_discrete_space_from_impl<DDim>(
typename NonUniformBSplines<typename DDim::continuous_dimension_type, DDim::degree()>::
template Impl<DDim, Kokkos::HostSpace>(breaks_begin, breaks_end));
}

/** @brief Constructs a non uniform spline basis using break points contained in a std::vector.
* The associated static attributes are available after this call.
*
* @tparam DDim The name of the discrete dimension
* @param breaks The std::vector of the coordinates of break points.
*/
template <class DDim>
void create_non_uniform_bsplines(
std::vector<ddc::Coordinate<typename DDim::continuous_dimension_type>> const& breaks)
{
create_non_uniform_bsplines<DDim>(breaks.begin(), breaks.end());
}

template <class DDim>
struct is_non_uniform_bsplines : public std::is_base_of<detail::NonUniformBSplinesBase, DDim>::type
{
Expand Down Expand Up @@ -439,7 +479,7 @@ NonUniformBSplines<CDim, D>::Impl<DDim, MemorySpace>::Impl(
knots[degree() + npoints() - 1 + i] = rmax;
}
}
ddc::init_discrete_space<knot_discrete_dimension_type>(knots);
ddc::create_non_uniform_point_sampling<knot_discrete_dimension_type>(knots);
}

template <class CDim, std::size_t D>
Expand Down
Loading
Loading