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.

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

C++17: if constexpr for templates in ShapeFactors #2659

Merged
Merged
Changes from all commits
Commits
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
234 changes: 84 additions & 150 deletions Source/Particles/ShapeFactors.H
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
/* Copyright 2019 Maxence Thevenet, Michael Rowan
/* Copyright 2019-2021 Maxence Thevenet, Michael Rowan, Luca Fedeli, Axel Huebl
*
* This file is part of WarpX.
*
Expand All @@ -7,10 +7,14 @@
#ifndef SHAPEFACTORS_H_
#define SHAPEFACTORS_H_

#include <AMReX.H>
#include <AMReX_GpuQualifiers.H>


/**
ax3l marked this conversation as resolved.
Show resolved Hide resolved
* Compute shape factor and return index of leftmost cell where
* particle writes.
* Specialized templates are defined below for orders 0 to 3.
* Specializations are defined for orders 0 to 3 (using "if constexpr").
* Shape factor functors may be evaluated with double arguments
* in current deposition to ensure that current deposited by
* particles that move only a small distance is still resolved.
Expand All @@ -23,170 +27,100 @@ struct Compute_shape_factor
{
template< typename T >
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
ax3l marked this conversation as resolved.
Show resolved Hide resolved
int operator()(T* const /*sx*/, T /*xint*/) const { return 0; }
};

/**
* Compute shape factor and return index of leftmost cell where
* particle writes.
* Specialization for order 0
*/
template <>
struct Compute_shape_factor< 0 >
{
template< typename T >
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
int operator()(T* const sx, T xmid) const
int operator()(
T* const sx,
T xmid) const
{
const auto j = static_cast<int>(xmid + T(0.5));
sx[0] = T(1.0);
return j;
if constexpr (depos_order == 0){
const auto j = static_cast<int>(xmid + T(0.5));
sx[0] = T(1.0);
return j;
}
else if constexpr (depos_order == 1){
const auto j = static_cast<int>(xmid);
const T xint = xmid - T(j);
sx[0] = T(1.0) - xint;
sx[1] = xint;
return j;
}
else if constexpr (depos_order == 2){
const auto j = static_cast<int>(xmid + T(0.5));
const T xint = xmid - T(j);
sx[0] = T(0.5)*(T(0.5) - xint)*(T(0.5) - xint);
sx[1] = T(0.75) - xint*xint;
sx[2] = T(0.5)*(T(0.5) + xint)*(T(0.5) + xint);
// index of the leftmost cell where particle deposits
return j-1;
}
else if constexpr (depos_order == 3){
const auto j = static_cast<int>(xmid);
const T xint = xmid - T(j);
sx[0] = (T(1.0))/(T(6.0))*(T(1.0) - xint)*(T(1.0) - xint)*(T(1.0) - xint);
sx[1] = (T(2.0))/(T(3.0)) - xint*xint*(T(1.0) - xint/(T(2.0)));
sx[2] = (T(2.0))/(T(3.0)) - (T(1.0) - xint)*(T(1.0) - xint)*(T(1.0) - T(0.5)*(T(1.0) - xint));
sx[3] = (T(1.0))/(T(6.0))*xint*xint*xint;
// index of the leftmost cell where particle deposits
return j-1;
}
else{
amrex::Abort("Unknown particle shape selected in Compute_shape_factor");
amrex::ignore_unused(sx, xmid);
ax3l marked this conversation as resolved.
Show resolved Hide resolved
return 0;
ax3l marked this conversation as resolved.
Show resolved Hide resolved
}
}
};

/**
* Compute shape factor and return index of leftmost cell where
* particle writes.
* Specialization for order 1
*/
template <>
struct Compute_shape_factor< 1 >
{
template< typename T >
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
int operator()(T* const sx, T xmid) const
{
const auto j = static_cast<int>(xmid);
const T xint = xmid - T(j);
sx[0] = T(1.0) - xint;
sx[1] = xint;
return j;
}
};

/**
* Compute shape factor and return index of leftmost cell where
* particle writes.
* Specialization for order 2
*/
template <>
struct Compute_shape_factor< 2 >
{
template< typename T >
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
int operator()(T* const sx, T xmid) const
{
const auto j = static_cast<int>(xmid + T(0.5));
const T xint = xmid - T(j);
sx[0] = T(0.5)*(T(0.5) - xint)*(T(0.5) - xint);
sx[1] = T(0.75) - xint*xint;
sx[2] = T(0.5)*(T(0.5) + xint)*(T(0.5) + xint);
// index of the leftmost cell where particle deposits
return j-1;
}
};

/**
* Compute shape factor and return index of leftmost cell where
* particle writes.
* Specialization for order 3
*/
template <>
struct Compute_shape_factor< 3 >
{
template< typename T >
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
int operator()(T* const sx, T xmid) const
{
const auto j = static_cast<int>(xmid);
const T xint = xmid - T(j);
sx[0] = (T(1.0))/(T(6.0))*(T(1.0) - xint)*(T(1.0) - xint)*(T(1.0) - xint);
sx[1] = (T(2.0))/(T(3.0)) - xint*xint*(T(1.0) - xint/(T(2.0)));
sx[2] = (T(2.0))/(T(3.0)) - (T(1.0) - xint)*(T(1.0) - xint)*(T(1.0) - T(0.5)*(T(1.0) - xint));
sx[3] = (T(1.0))/(T(6.0))*xint*xint*xint;
// index of the leftmost cell where particle deposits
return j-1;
}
};

/**
* Compute shifted shape factor and return index of leftmost cell where
* particle writes, for Esirkepov algorithm.
* Specialized templates are defined below for orders 1, 2 and 3.
* Specializations are defined below for orders 1, 2 and 3 (using "if constexpr").
*/
template <int depos_order>
struct Compute_shifted_shape_factor
{
template< typename T >
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
int operator()(T* const sx, const T x_old, const int i_new) const;
};

/**
* Compute shifted shape factor and return index of leftmost cell where
* particle writes, for Esirkepov algorithm.
* Specialization for order 1
*/
template <>
struct Compute_shifted_shape_factor< 1 >
{
template< typename T >
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
int operator()(T* const sx, const T x_old, const int i_new) const
{
const auto i = static_cast<int>(x_old);
const int i_shift = i - i_new;
const T xint = x_old - T(i);
sx[1+i_shift] = T(1.0) - xint;
sx[2+i_shift] = xint;
return i;
}
};

/**
* Compute shifted shape factor and return index of leftmost cell where
* particle writes, for Esirkepov algorithm.
* Specialization for order 2
*/
template <>
struct Compute_shifted_shape_factor< 2 >
{
template< typename T >
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
int operator()(T* const sx, const T x_old, const int i_new) const
{
const auto i = static_cast<int>(x_old + T(0.5));
const int i_shift = i - (i_new + 1);
const T xint = x_old - T(i);
sx[1+i_shift] = T(0.5)*(T(0.5) - xint)*(T(0.5) - xint);
sx[2+i_shift] = T(0.75) - xint*xint;
sx[3+i_shift] = T(0.5)*(T(0.5) + xint)*(T(0.5) + xint);
// index of the leftmost cell where particle deposits
return i - 1;
}
};

/**
* Compute shifted shape factor and return index of leftmost cell where
* particle writes, for Esirkepov algorithm.
* Specialization for order 3
*/
template <>
struct Compute_shifted_shape_factor< 3 >
{
template< typename T >
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
int operator()(T* const sx, const T x_old, const int i_new) const
int operator()(
T* const sx,
const T x_old,
const int i_new) const
{
const auto i = static_cast<int>(x_old);
const int i_shift = i - (i_new + 1);
const T xint = x_old - i;
sx[1+i_shift] = (T(1.0))/(T(6.0))*(T(1.0) - xint)*(T(1.0) - xint)*(T(1.0) - xint);
sx[2+i_shift] = (T(2.0))/(T(3.0)) - xint*xint*(T(1.0) - xint/(T(2.0)));
sx[3+i_shift] = (T(2.0))/(T(3.0)) - (T(1.0) - xint)*(T(1.0) - xint)*(T(1.0) - T(0.5)*(T(1.0) - xint));
sx[4+i_shift] = (T(1.0))/(T(6.0))*xint*xint*xint;
// index of the leftmost cell where particle deposits
return i - 1;
if constexpr (depos_order == 1){
const auto i = static_cast<int>(x_old);
const int i_shift = i - i_new;
const T xint = x_old - T(i);
sx[1+i_shift] = T(1.0) - xint;
sx[2+i_shift] = xint;
return i;
}
else if constexpr (depos_order == 2){
const auto i = static_cast<int>(x_old + T(0.5));
const int i_shift = i - (i_new + 1);
const T xint = x_old - T(i);
sx[1+i_shift] = T(0.5)*(T(0.5) - xint)*(T(0.5) - xint);
sx[2+i_shift] = T(0.75) - xint*xint;
sx[3+i_shift] = T(0.5)*(T(0.5) + xint)*(T(0.5) + xint);
// index of the leftmost cell where particle deposits
return i - 1;
}
else if constexpr (depos_order == 3){
const auto i = static_cast<int>(x_old);
const int i_shift = i - (i_new + 1);
const T xint = x_old - i;
sx[1+i_shift] = (T(1.0))/(T(6.0))*(T(1.0) - xint)*(T(1.0) - xint)*(T(1.0) - xint);
sx[2+i_shift] = (T(2.0))/(T(3.0)) - xint*xint*(T(1.0) - xint/(T(2.0)));
sx[3+i_shift] = (T(2.0))/(T(3.0)) - (T(1.0) - xint)*(T(1.0) - xint)*(T(1.0) - T(0.5)*(T(1.0) - xint));
sx[4+i_shift] = (T(1.0))/(T(6.0))*xint*xint*xint;
// index of the leftmost cell where particle deposits
return i - 1;
}
else{
amrex::Abort("Unknown particle shape selected in Compute_shifted_shape_factor");
amrex::ignore_unused(sx, x_old, i_new);
return 0;
ax3l marked this conversation as resolved.
Show resolved Hide resolved
}
}
};

Expand Down