AMReX_IntVect.H

Go to the documentation of this file.

#ifndef AMREX_INTVECT_H_
#define AMREX_INTVECT_H_
#include <AMReX_Config.H>
 
#include <AMReX_INT.H>
#include <AMReX_SPACE.H>
#include <AMReX_Array.H>
#include <AMReX_Vector.H>
#include <AMReX_Dim3.H>
#include <AMReX_BLassert.H>
#include <AMReX_Extension.H>
#include <AMReX_GpuQualifiers.H>
#include <AMReX_Math.H>
#include <AMReX_Tuple.H>
#include <AMReX_TypeTraits.H>
 
#include <iosfwd>
#include <cstdlib>
#include <cmath>
#include <limits>
#include <climits>
#include <algorithm>
 
namespace amrex {
 
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
int coarsen (int i, int ratio) noexcept
{
    switch (ratio) {
    case  1: return i;
    case  2: return (i<0) ? -std::abs(i+1)/2    -1 : i/2;
    case  4: return (i<0) ? -std::abs(i+1)/4    -1 : i/4;
    default: return (i<0) ? -std::abs(i+1)/ratio-1 : i/ratio;
    }
}
 
template<int dim>
class IntVectND
{
public:
    static_assert(dim >= 1, "The number of dimensions of IntVectND must be positive");
 
    struct shift_hasher {
        std::size_t operator()(const IntVectND<dim>& vec) const noexcept
        {
            static constexpr unsigned shift1 = sizeof(size_t)>=8 ? 20 : 10;
            static constexpr unsigned shift2 = sizeof(size_t)>=8 ? 40 : 20;
            if constexpr (dim == 1) {
                amrex::ignore_unused(shift1);
                amrex::ignore_unused(shift2);
                return static_cast<std::size_t>(vec[0]);
            } else if constexpr (dim == 2) {
                amrex::ignore_unused(shift2);
                return static_cast<std::size_t>(vec[0]) ^
                      (static_cast<std::size_t>(vec[1]) << shift1);
            } else if constexpr (dim == 3) {
                return static_cast<std::size_t>(vec[0]) ^
                      (static_cast<std::size_t>(vec[1]) << shift1) ^
                      (static_cast<std::size_t>(vec[2]) << shift2);
            } else {
                std::size_t seed = dim;
                // hash function from
                // https://stackoverflow.com/questions/20511347/a-good-hash-function-for-a-vector
                for (int i=0; i<dim; ++i) {
                    auto x = static_cast<unsigned int>(vec[i]);
                    x = ((x >> 16) ^ x) * 0x45d9f3b;
                    x = ((x >> 16) ^ x) * 0x45d9f3b;
                    x = (x >> 16) ^ x;
                    seed ^= x + 0x9e3779b9 + (seed << 6) + (seed >> 2);
                }
                return seed;
            }
        }
    };
 
    using hasher = shift_hasher;
 
 
    constexpr IntVectND () noexcept {} // cannot use = default due to Clang bug // NOLINT
 
    template <class...Args,
        std::enable_if_t<
            (sizeof...(Args)+2 == dim) &&
            IsConvertible_v<int, Args...>,
        int> = 0>
    AMREX_GPU_HOST_DEVICE
    constexpr IntVectND (int i, int j, Args...ks) noexcept : vect{i, j, static_cast<int>(ks)...} {}
 
    AMREX_GPU_HOST_DEVICE
    explicit constexpr IntVectND (int s) noexcept {
        for (int i=0; i<dim; ++i) {
            vect[i] = s;
        }
    }
 
    AMREX_GPU_HOST_DEVICE
    explicit IntVectND (const int* a) noexcept {
        for (int i=0; i<dim; ++i) {
            vect[i] = a[i];
        }
    }
 
    explicit IntVectND (const Vector<int>& a) noexcept {
        BL_ASSERT(a.size() == dim);
        for (int i=0; i<dim; ++i) {
            vect[i] = a[i];
        }
    }
 
    explicit IntVectND (const Array<int,dim>& a) noexcept {
        for (int i=0; i<dim; ++i) {
            vect[i] = a[i];
        }
    }
 
    template <int N=dim, std::enable_if_t<( 1<=N && N<=3 ), int> = 0>
    explicit constexpr IntVectND (Dim3 const& a) noexcept {
        vect[0] = a.x;
        if constexpr (dim >= 2) {
            vect[1] = a.y;
        }
        if constexpr (dim == 3) {
            vect[2] = a.z;
        }
    }
 
    // dtor, copy-ctor, copy-op=, move-ctor, and move-op= are compiler generated.
 
    template <int N=dim, std::enable_if_t<( 1<=N && N<=3 ), int> = 0>
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    Dim3 dim3 () const noexcept {
        if constexpr (dim == 1) {
            return Dim3{vect[0],0,0};
        } else if constexpr (dim == 2) {
            return Dim3{vect[0],vect[1],0};
        } else {
            return Dim3{vect[0],vect[1],vect[2]};
        }
    }
 
    template <int N=dim, std::enable_if_t<( 1<=N && N<=3 ), int> = 0>
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    Dim3 dim3 ([[maybe_unused]] int fill_extra) const noexcept {
        if constexpr (dim == 1) {
            return Dim3{vect[0],fill_extra,fill_extra};
        } else if constexpr (dim == 2) {
            return Dim3{vect[0],vect[1],fill_extra};
        } else {
            return Dim3{vect[0],vect[1],vect[2]};
        }
    }
 
#if __cplusplus >= 201402L
    template< typename T = int >
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    Array<T, dim>
    toArray () const noexcept {
        Array<T, dim> ret {};
        for (int i=0; i<dim; ++i) {
            ret[i] = T(vect[i]);
        }
        return ret;
    }
#endif
 
 
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    int sum () const noexcept
    {
        int retval = vect[0];
        for (int i=1; i<dim; ++i) {
            retval += vect[i];
        }
        return retval;
    }
 
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    int max () const noexcept
    {
        int retval = vect[0];
        for (int i=1; i<dim; ++i) {
            retval = retval > vect[i] ? retval : vect[i];
        }
        return retval;
    }
 
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    int min () const noexcept
    {
        int retval = vect[0];
        for (int i=1; i<dim; ++i) {
            retval = retval < vect[i] ? retval : vect[i];
        }
        return retval;
    }
 
    //return coordinate with largest value
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    int maxDir(bool a_doAbsValue) const noexcept;
 
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE constexpr
    int& operator[] (int i) noexcept { BL_ASSERT(i>=0 && i < dim); return vect[i]; }
 
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE constexpr
    const int& operator[] (int i) const noexcept { BL_ASSERT(i>=0 && i < dim); return vect[i]; }
 
    template<std::size_t i>
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE constexpr
    int& get () noexcept {static_assert(0<=i && i<dim); return vect[i];}
 
    template<std::size_t i>
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE constexpr
    const int& get () const noexcept {static_assert(0<=i && i<dim); return vect[i];}
 
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE constexpr
    int* begin () noexcept { return &vect[0]; }
 
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE constexpr
    const int* begin () const noexcept { return &vect[0]; }
 
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE constexpr
    int* end () noexcept { return &vect[dim]; }
 
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE constexpr
    const int* end () const noexcept { return &vect[dim]; }
 
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND& setVal (int i, int val) noexcept
    {
        BL_ASSERT(i>=0 && i<dim); vect[i] = val; return *this;
    }
 
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    const int* getVect () const& noexcept { return vect; }
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    int*       getVect () &      noexcept { return vect; }
    AMREX_GPU_HOST_DEVICE
    int*       getVect () && = delete;
 
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    bool operator== (int val) const noexcept
    {
        bool retval = vect[0] == val;
        for (int i=1; i<dim; ++i) {
            retval = retval && vect[i] == val;
        }
        return retval;
    }
 
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    bool operator!= (int val) const noexcept
    {
        bool retval = vect[0] != val;
        for (int i=1; i<dim; ++i) {
            retval = retval || vect[i] != val;
        }
        return retval;
    }
 
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    bool operator== (const IntVectND<dim>& rhs) const noexcept
    {
        bool retval = vect[0] == rhs[0];
        for (int i=1; i<dim; ++i) {
            retval = retval && vect[i] == rhs[i];
        }
        return retval;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    bool operator!= (const IntVectND<dim>& rhs) const noexcept
    {
        bool retval = vect[0] != rhs[0];
        for (int i=1; i<dim; ++i) {
            retval = retval || vect[i] != rhs[i];
        }
        return retval;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    bool operator< (const IntVectND<dim>& rhs) const noexcept
    {
        for (int i=dim-1; i>=0; --i) {
            if (vect[i] < rhs[i]) {
                return true;
            } else if (vect[i] > rhs[i]) {
                return false;
            }
        }
        return false;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    bool operator<= (const IntVectND<dim>& rhs) const noexcept
    {
        return !(rhs < *this);
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    bool operator> (const IntVectND<dim>& rhs) const noexcept
    {
        return rhs < *this;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    bool operator>= (const IntVectND<dim>& rhs) const noexcept
    {
        return !(*this < rhs);
    }
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    bool allLT (const IntVectND<dim>& rhs) const noexcept
    {
        bool retval = vect[0] < rhs[0];
        for (int i=1; i<dim; ++i) {
            retval = retval && vect[i] < rhs[i];
        }
        return retval;
    }
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    bool allLT (int rhs) const noexcept
    {
        bool retval = vect[0] < rhs;
        for (int i=1; i<dim; ++i) {
            retval = retval && vect[i] < rhs;
        }
        return retval;
    }
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    bool allLE (const IntVectND<dim>& rhs) const noexcept
    {
        bool retval = vect[0] <= rhs[0];
        for (int i=1; i<dim; ++i) {
            retval = retval && vect[i] <= rhs[i];
        }
        return retval;
    }
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    bool allLE (int rhs) const noexcept
    {
        bool retval = vect[0] <= rhs;
        for (int i=1; i<dim; ++i) {
            retval = retval && vect[i] <= rhs;
        }
        return retval;
    }
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    bool allGT (const IntVectND<dim>& rhs) const noexcept
    {
        bool retval = vect[0] > rhs[0];
        for (int i=1; i<dim; ++i) {
            retval = retval && vect[i] > rhs[i];
        }
        return retval;
    }
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    bool allGT (int rhs) const noexcept
    {
        bool retval = vect[0] > rhs;
        for (int i=1; i<dim; ++i) {
            retval = retval && vect[i] > rhs;
        }
        return retval;
    }
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    bool allGE (const IntVectND<dim>& rhs) const noexcept
    {
        bool retval = vect[0] >= rhs[0];
        for (int i=1; i<dim; ++i) {
            retval = retval && vect[i] >= rhs[i];
        }
        return retval;
    }
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    bool allGE (int rhs) const noexcept
    {
        bool retval = vect[0] >= rhs;
        for (int i=1; i<dim; ++i) {
            retval = retval && vect[i] >= rhs;
        }
        return retval;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim> operator+ () const noexcept { return *this; }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim> operator- () const noexcept {
        IntVectND<dim> retval(0);
        for (int i=0; i<dim; ++i) {
            retval[i] = -vect[i];
        }
        return retval;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim>& operator+= (int s) noexcept
    {
        for (int i=0; i<dim; ++i) {
            vect[i] += s;
        }
        return *this;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim>& operator+= (const IntVectND<dim>& p) noexcept
    {
        for (int i=0; i<dim; ++i) {
            vect[i] += p[i];
        }
        return *this;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim>& operator*= (int s) noexcept
    {
        for (int i=0; i<dim; ++i) {
            vect[i] *= s;
        }
        return *this;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim>& operator*= (const IntVectND<dim>& p) noexcept
    {
        for (int i=0; i<dim; ++i) {
            vect[i] *= p[i];
        }
        return *this;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim>& operator/= (int s) noexcept
    {
        for (int i=0; i<dim; ++i) {
            vect[i] /= s;
        }
        return *this;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim>& operator/= (const IntVectND<dim>& p) noexcept
    {
        for (int i=0; i<dim; ++i) {
            vect[i] /= p[i];
        }
        return *this;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim>& operator-= (int s) noexcept
    {
        for (int i=0; i<dim; ++i) {
            vect[i] -= s;
        }
        return *this;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim>& operator-= (const IntVectND<dim>& p) noexcept
    {
        for (int i=0; i<dim; ++i) {
            vect[i] -= p[i];
        }
        return *this;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim> operator+ (const IntVectND<dim>& p) const noexcept
    {
        IntVectND<dim> retval = *this;
        return retval += p;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim> operator+ (int s) const noexcept
    {
        IntVectND<dim> retval = *this;
        return retval += s;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim> operator- (const IntVectND<dim>& p) const noexcept
    {
        IntVectND<dim> retval = *this;
        return retval -= p;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim> operator- (int s) const noexcept
    {
        IntVectND<dim> retval = *this;
        return retval -= s;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim> operator* (const IntVectND<dim>& p) const noexcept
    {
        IntVectND<dim> retval = *this;
        return retval *= p;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim> operator* (int s) const noexcept
    {
        IntVectND<dim> retval = *this;
        return retval *= s;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim> operator/ (const IntVectND<dim>& p) const noexcept
    {
        IntVectND<dim> retval = *this;
        return retval /= p;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim> operator/ (int s) const noexcept
    {
        IntVectND<dim> retval = *this;
        return retval /= s;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim>& min (const IntVectND<dim>& p) noexcept
    {
        for (int i=0; i<dim; ++i) {
            vect[i] = (vect[i] < p.vect[i] ? vect[i] : p.vect[i]);
        }
        return *this;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim>& max (const IntVectND<dim>& p) noexcept
    {
        for (int i=0; i<dim; ++i) {
            vect[i] = (vect[i] > p.vect[i] ? vect[i] : p.vect[i]);
        }
        return *this;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim>& scale (int s) noexcept {
        for (int i=0; i<dim; ++i) {
            vect[i] *= s;
        }
        return *this;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim>& reflect (int ref_ix, int idir) noexcept
    {
        BL_ASSERT(idir >= 0 && idir < dim);
        vect[idir] = -vect[idir] + 2*ref_ix;
        return *this;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim>& shift (int coord, int s) noexcept
    {
        BL_ASSERT(coord >= 0 && coord < dim); vect[coord] += s; return *this;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim>& shift (const IntVectND<dim>& iv) noexcept { *this += iv; return *this; }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim>& diagShift (int s) noexcept
    {
        for (int i=0; i<dim; ++i) {
            vect[i] += s;
        }
        return *this;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim>& coarsen (const IntVectND<dim>& p) noexcept;
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<dim>& coarsen (int p) noexcept;
 
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    static constexpr IntVectND<dim> TheZeroVector () noexcept {
        return IntVectND<dim>(0);
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    static constexpr IntVectND<dim> TheUnitVector () noexcept {
        return IntVectND<dim>(1);
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    static constexpr IntVectND<dim> TheDimensionVector (int d) noexcept {
        IntVectND<dim> retval(0);
        retval[d] = 1;
        return retval;
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    static constexpr IntVectND<dim> TheNodeVector () noexcept {
        return IntVectND<dim>(1);
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    static constexpr IntVectND<dim> TheCellVector () noexcept {
        return IntVectND<dim>(0);
    }
 
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    static constexpr IntVectND<dim> TheMaxVector () noexcept {
        return IntVectND<dim>(std::numeric_limits<int>::max());
    }
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    static constexpr IntVectND<dim> TheMinVector () noexcept {
        return IntVectND<dim>(std::numeric_limits<int>::lowest());
    }
 
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    static constexpr std::size_t size () noexcept {
        return static_cast<std::size_t>(dim);
    }
 
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    static constexpr int isize () noexcept {
        return dim;
    }
 
    using value_type = int;
 
    template<int new_dim>
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<new_dim> shrink () const noexcept {
        static_assert(new_dim <= dim);
        return IntVectND<new_dim>(this->begin());
    }
 
    template<int new_dim>
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<new_dim> expand (int fill_extra=0) const noexcept {
        static_assert(new_dim >= dim);
        IntVectND<new_dim> retval(fill_extra);
        for (int i=0; i<dim; ++i) {
            retval[i] = vect[i];
        }
        return retval;
    }
 
    template<int new_dim>
    [[nodiscard]] AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
    IntVectND<new_dim> resize (int fill_extra=0) const noexcept {
        if constexpr (new_dim > dim) {
            return expand<new_dim>(fill_extra);
        } else {
            return shrink<new_dim>();
        }
    }
 
    static const IntVectND<dim> Zero;
 
    static const IntVectND<dim> Unit;
 
private:
 
    int vect[dim] = {};
};
 
template <int dim>
inline constexpr const IntVectND<dim> IntVectND<dim>::Zero{0};
 
template <int dim>
inline constexpr const IntVectND<dim> IntVectND<dim>::Unit{1};
 
// Template deduction guide for IntVectND
template<std::size_t dim>
AMREX_GPU_HOST_DEVICE // __device__ for HIP
IntVectND(const Array<int, dim>&) -> IntVectND<dim>;
 
// Template deduction guide for IntVectND
template <class...Args,
    std::enable_if_t<
        IsConvertible_v<int, Args...>,
    int> = 0>
AMREX_GPU_HOST_DEVICE // __device__ for HIP
IntVectND(int, int, Args...) -> IntVectND<sizeof...(Args)+2>;
 
using IntVect = IntVectND<AMREX_SPACEDIM>;
 
template<int dim>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
IntVectND<dim>&
IntVectND<dim>::coarsen (int s) noexcept
{
    BL_ASSERT(s > 0);
    switch (s) {
    case 1:
        break;
    case 2:
        for (int i=0; i<dim; ++i) {
            vect[i] = (vect[i]<0) ? -std::abs(vect[i]+1)/2-1 : vect[i]/2;
        }
        break;
    case 4:
        for (int i=0; i<dim; ++i) {
            vect[i] = (vect[i]<0) ? -std::abs(vect[i]+1)/4-1 : vect[i]/4;
        }
        break;
    default:
        for (int i=0; i<dim; ++i) {
            vect[i] = (vect[i]<0) ? -std::abs(vect[i]+1)/s-1 : vect[i]/s;
        }
    }
    return *this;
}
 
template<int dim>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
IntVectND<dim>&
IntVectND<dim>::coarsen (const IntVectND<dim>& p) noexcept
{
    BL_ASSERT(p.allGT(0));
    for (int i=0; i<dim; ++i) {
        vect[i] = amrex::coarsen(vect[i], p.vect[i]);
    }
    return *this;
}
 
template<int dim>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
int
IntVectND<dim>::maxDir(bool a_doAbsValue) const noexcept
{
    int retval = 0;
    if(a_doAbsValue)
    {
        int maxval = std::abs((*this)[0]);
        for(int idir = 1; idir < dim; idir++)
        {
            int curval = std::abs((*this)[idir]);
            if(curval > maxval)
            {
                maxval = curval;
                retval = idir;
            }
        }
    }
    else
    {
        int maxval = (*this)[0];
        for(int idir = 1; idir < dim; idir++)
        {
            int curval = (*this)[idir];
            if(curval > maxval)
            {
                maxval = curval;
                retval = idir;
            }
        }
    }
    return retval;
}
 
template<int dim>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
IntVectND<dim> operator+ (int s, const IntVectND<dim>& p) noexcept
{
    IntVectND<dim> retval = p;
    for (int i=0; i<dim; ++i) {
        retval[i] = s + retval[i];
    }
    return retval;
}
 
template<int dim>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
AMREX_GPU_HOST_DEVICE
IntVectND<dim> operator- (int s, const IntVectND<dim>& p) noexcept
{
    IntVectND<dim> retval = p;
    for (int i=0; i<dim; ++i) {
        retval[i] = s - retval[i];
    }
    return retval;
}
 
template<int dim>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
IntVectND<dim> operator* (int s, const IntVectND<dim>& p) noexcept
{
    IntVectND<dim> retval = p;
    for (int i=0; i<dim; ++i) {
        retval[i] = s * retval[i];
    }
    return retval;
}
 
template<int dim>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
IntVectND<dim>
min (const IntVectND<dim>& p1, const IntVectND<dim>& p2) noexcept
{
    IntVectND<dim> p(p1);
    p.min(p2);
    return p;
}
 
template<int dim>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
IntVectND<dim>
elemwiseMin (const IntVectND<dim>& p1, const IntVectND<dim>& p2) noexcept
{
    IntVectND<dim> p(p1);
    p.min(p2);
    return p;
}
 
template<int dim>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
IntVectND<dim>
max (const IntVectND<dim>& p1, const IntVectND<dim>& p2) noexcept
{
    IntVectND<dim> p(p1);
    p.max(p2);
    return p;
}
 
template<int dim>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
IntVectND<dim>
elemwiseMax (const IntVectND<dim>& p1, const IntVectND<dim>& p2) noexcept
{
    IntVectND<dim> p(p1);
    p.max(p2);
    return p;
}
 
template<int dim = AMREX_SPACEDIM>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
IntVectND<dim>
BASISV (int dir) noexcept
{
    BL_ASSERT(dir >= 0 && dir < dim);
    IntVectND<dim> tmp(0);
    tmp[dir] = 1;
    return tmp;
}
 
template<int dim>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
IntVectND<dim>
scale (const IntVectND<dim>& p, int s) noexcept
{
    return IntVectND<dim>(p).scale(s);
}
 
template<int dim>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
IntVectND<dim>
reflect (const IntVectND<dim>& a, int ref_ix, int idir) noexcept
{
    return IntVectND<dim>(a).reflect(ref_ix, idir);
}
 
template<int dim>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
IntVectND<dim>
diagShift (const IntVectND<dim>& p, int s) noexcept
{
    return p + s;
}
 
template<int dim>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
IntVectND<dim>
coarsen (const IntVectND<dim>& p, int s) noexcept
{
    BL_ASSERT(s > 0);
    IntVectND<dim> v = p;
    v.coarsen(s);
    return v;
}
 
template<int dim>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
IntVectND<dim>
coarsen (const IntVectND<dim>& p1, const IntVectND<dim>& p2) noexcept
{
    IntVectND<dim> v = p1;
    v.coarsen(p2);
    return v;
}
 
template<int dim, std::enable_if_t<( 1<=dim && dim<=3 ), int> = 0>
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
Dim3 refine (Dim3 const& coarse, IntVectND<dim> const& ratio) noexcept
{
    if constexpr (dim == 1) {
        return Dim3{coarse.x*ratio[0], coarse.y, coarse.z};
    } else if constexpr (dim == 2) {
        return Dim3{coarse.x*ratio[0], coarse.y*ratio[1], coarse.z};
    } else {
        return Dim3{coarse.x*ratio[0], coarse.y*ratio[1], coarse.z*ratio[2]};
    }
}
 
template<int dim, std::enable_if_t<( 1<=dim && dim<=3 ), int> = 0>
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
Dim3 coarsen (Dim3 const& fine, IntVectND<dim> const& ratio) noexcept
{
    if constexpr (dim == 1) {
        return Dim3{amrex::coarsen(fine.x, ratio[0]),
                    fine.y,
                    fine.z};
    } else if constexpr (dim == 2) {
        return Dim3{amrex::coarsen(fine.x, ratio[0]),
                    amrex::coarsen(fine.y, ratio[1]),
                    fine.z};
    } else {
        return Dim3{amrex::coarsen(fine.x, ratio[0]),
                    amrex::coarsen(fine.y, ratio[1]),
                    amrex::coarsen(fine.z, ratio[2])};
    }
}
 
namespace detail {
    std::ostream& int_vector_write (std::ostream& os, const int* iv, int dim);
    std::istream& int_vector_read (std::istream& is, int* iv, int dim);
 
    template<int dim>
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE constexpr
    void IntVectCat_imp (int*& dst, const IntVectND<dim>& src) noexcept {
        for (int i=0; i<dim; ++i) {
            dst[i] = src[i];
        }
        dst += dim;
    }
 
    template<int dim>
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE constexpr
    void IntVectSplit_imp2 (IntVectND<dim>& dst, const int*& src) noexcept {
        for (int i=0; i<dim; ++i) {
            dst[i] = src[i];
        }
        src += dim;
    }
 
    template<class T, std::size_t...Ns>
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE constexpr
    T IntVectSplit_imp (T& retval, std::index_sequence<Ns...>, const int * src) noexcept {
        (IntVectSplit_imp2(amrex::get<Ns>(retval), src), ...);
        return retval;
    }
 
    template<int...dims>
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE constexpr
    int get_sum () {
        return (0 + ... + dims);
    }
}
 
template<int dim>
std::ostream&
operator<< (std::ostream& os, const IntVectND<dim>& iv)
{
    return detail::int_vector_write(os, iv.begin(), dim);
}
 
template<int dim>
std::istream&
operator>> (std::istream& is, IntVectND<dim>& iv)
{
    return detail::int_vector_read(is, iv.begin(), dim);
}
 
template<int d, int...dims>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
constexpr IntVectND<detail::get_sum<d, dims...>()>
IntVectCat (const IntVectND<d>& v, const IntVectND<dims>&...vects) noexcept {
    IntVectND<detail::get_sum<d, dims...>()> retval (0);
    int* dst = retval.begin();
    detail::IntVectCat_imp(dst, v);
    (detail::IntVectCat_imp(dst, vects), ...);
    return retval;
}
 
template<int d, int...dims>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
constexpr GpuTuple<IntVectND<d>, IntVectND<dims>...>
IntVectSplit (const IntVectND<detail::get_sum<d, dims...>()>& v) noexcept {
    GpuTuple<IntVectND<d>, IntVectND<dims>...> retval;
    return detail::IntVectSplit_imp(retval,
                                    std::make_index_sequence<1 + sizeof...(dims)>(),
                                    v.begin());
}
 
template<int new_dim, int old_dim>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
constexpr IntVectND<new_dim>
IntVectShrink (const IntVectND<old_dim>& iv) noexcept {
    return iv.template shrink<new_dim>();
}
 
template<int new_dim, int old_dim>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
constexpr IntVectND<new_dim>
IntVectExpand (const IntVectND<old_dim>& iv, int fill_extra=0) noexcept {
    return iv.template expand<new_dim>(fill_extra);
}
 
template<int new_dim, int old_dim>
AMREX_GPU_HOST_DEVICE
AMREX_FORCE_INLINE
constexpr IntVectND<new_dim>
IntVectResize (const IntVectND<old_dim>& iv, int fill_extra=0) noexcept {
    return iv.template resize<new_dim>(fill_extra);
}
 
} // namespace amrex
 
// Spcialize std::tuple_size for IntVectND. Used by structured bindings.
template<int dim>
struct std::tuple_size<amrex::IntVectND<dim>> {
    static constexpr std::size_t value = dim;
};
 
// Spcialize std::tuple_element for IntVectND. Used by structured bindings.
template<std::size_t s, int dim>
struct std::tuple_element<s, amrex::IntVectND<dim>> {
    using type = int;
};
 
#endif /*AMREX_INTVECT_H*/