single-array-storage-base.hpp

/*
 * single-array-storage-base.hpp
 *
 * Created on: 07 Nov 2022
 * Author: Holger Schmitz
 * Email: holger@notjustphysics.com
 *
 * Copyright 2012-2022 Holger Schmitz
 *
 * This file is part of Schnek.
 *
 * Schnek is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 *
 * Schnek is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with Schnek.  If not, see <http://www.gnu.org/licenses/>.
 *
 */

#ifndef SCHNEK_GRID_GRIDSTORAGE_SINGLESTORAGEBASE_HPP_
#define SCHNEK_GRID_GRIDSTORAGE_SINGLESTORAGEBASE_HPP_

#include "../array.hpp"

namespace schnek
{
    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    class SingleArrayGridStorageBase : public AllocationPolicy<T, rank>
    {
    public:
        typedef Array<int, rank> IndexType;
        
        typedef Range<int, rank> RangeType;
    public:
        SingleArrayGridStorageBase();

        SingleArrayGridStorageBase(const SingleArrayGridStorageBase&) = default;

        SingleArrayGridStorageBase<T, rank, AllocationPolicy> &operator=(
            const SingleArrayGridStorageBase<T, rank, AllocationPolicy>&
        ) = default;

        T *getRawData() const { return this->data->ptr; }

        SCHNEK_INLINE const IndexType &getLo() const { return this->range.getLo(); }

        SCHNEK_INLINE const IndexType &getHi() const { return this->range.getHi(); }

        SCHNEK_INLINE const RangeType &getRange() const { return this->range; }

        SCHNEK_INLINE const IndexType &getDims() const { return this->dims; }

        SCHNEK_INLINE int getLo(int k) const { return this->range.getLo(k); }

        SCHNEK_INLINE int getHi(int k) const { return this->range.getHi(k); }

        SCHNEK_INLINE int getDims(int k) const { return this->dims[k]; }

        SCHNEK_INLINE int getSize() const { return this->size; }

        typedef T* storage_iterator;
        typedef const T* const_storage_iterator;

        SCHNEK_INLINE storage_iterator begin() { return this->data->ptr; }
        SCHNEK_INLINE storage_iterator end() { return this->data->ptr + this->size; }

        SCHNEK_INLINE const_storage_iterator cbegin() const { return this->data->ptr; }
        SCHNEK_INLINE const_storage_iterator cend() const { return this->data->ptr + this->size; }
    };

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    class SingleArrayGridCOrderStorageBase : public SingleArrayGridStorageBase<T, rank, AllocationPolicy>
    {
    private:
        T *data_fast;
    public:
        typedef SingleArrayGridStorageBase<T, rank, AllocationPolicy> BaseType;

        typedef typename BaseType::IndexType IndexType;

        typedef typename BaseType::RangeType RangeType;

        SingleArrayGridCOrderStorageBase();

        SingleArrayGridCOrderStorageBase(const SingleArrayGridCOrderStorageBase&);

        SingleArrayGridCOrderStorageBase(const IndexType &lo, const IndexType &hi);

        SingleArrayGridCOrderStorageBase(const RangeType &range);

        SingleArrayGridCOrderStorageBase<T, rank, AllocationPolicy> &operator=(
            const SingleArrayGridCOrderStorageBase<T, rank, AllocationPolicy> &
        ) = default;

        SCHNEK_INLINE T &get(const IndexType &index);

        SCHNEK_INLINE const T &get(const IndexType &index) const;

        void resize(const IndexType &low, const IndexType &high);

        void resize(const RangeType range);

        ptrdiff_t stride(size_t dim) const;
    private:
        void updateDataFast();

    };

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    class SingleArrayGridFortranOrderStorageBase : public SingleArrayGridStorageBase<T, rank, AllocationPolicy>
    {
    private:
        T *data_fast;
    public:
        typedef SingleArrayGridStorageBase<T, rank, AllocationPolicy> BaseType;

        typedef typename BaseType::IndexType IndexType;

        typedef typename BaseType::RangeType RangeType;

        SingleArrayGridFortranOrderStorageBase();

        SingleArrayGridFortranOrderStorageBase(const SingleArrayGridFortranOrderStorageBase&);
        
        SingleArrayGridFortranOrderStorageBase(const IndexType &lo, const IndexType &hi);

        SingleArrayGridFortranOrderStorageBase(const RangeType &range);

        SingleArrayGridFortranOrderStorageBase<T, rank, AllocationPolicy> &operator=(
            const SingleArrayGridFortranOrderStorageBase<T, rank, AllocationPolicy> &
        ) = default;

        SCHNEK_INLINE T &get(const IndexType &index);

        SCHNEK_INLINE const T &get(const IndexType &index) const;

        void resize(const IndexType &low, const IndexType &high);

        void resize(const RangeType range);

        ptrdiff_t stride(size_t dim) const;  
    private:
        void updateDataFast();  
    };

    //=================================================================
    //================== SingleArrayGridStorageBase ===================
    //=================================================================

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    SingleArrayGridStorageBase<T, rank, AllocationPolicy>::SingleArrayGridStorageBase()
        : AllocationPolicy<T, rank>()
    {
    }

    //=================================================================
    //=============== SingleArrayGridCOrderStorageBase ================
    //=================================================================

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    SingleArrayGridCOrderStorageBase<T, rank, AllocationPolicy>::SingleArrayGridCOrderStorageBase() 
        : BaseType(), data_fast(NULL) 
    {
        this->onUpdate([this](){ updateDataFast(); });
    }

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    SingleArrayGridCOrderStorageBase<T, rank, AllocationPolicy>::SingleArrayGridCOrderStorageBase(const SingleArrayGridCOrderStorageBase &other)
        : BaseType(other), data_fast(other.data_fast)
    {
        this->onUpdate([this](){ updateDataFast(); });
    }

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    SingleArrayGridCOrderStorageBase<T, rank, AllocationPolicy>::SingleArrayGridCOrderStorageBase(
        const IndexType &lo, 
        const IndexType &hi
    ) : BaseType(), data_fast(NULL)
    {
        this->onUpdate([this](){ updateDataFast(); });
        resize(lo, hi);
    }

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    SingleArrayGridCOrderStorageBase<T, rank, AllocationPolicy>::SingleArrayGridCOrderStorageBase(
        const RangeType &range
    ) : BaseType(), data_fast(NULL)
    {
        this->onUpdate([this](){ updateDataFast(); });
        resize(range.getLo(), range.getHi());
    }

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    SCHNEK_INLINE T &SingleArrayGridCOrderStorageBase<T, rank, AllocationPolicy>::get(const IndexType &index)
    {
        size_t pos = index[0];
        for (size_t i = 1; i < rank; ++i)
        {
            pos = index[i] + this->dims[i] * pos;
        }
        return this->data_fast[pos];
    }

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    SCHNEK_INLINE const T &SingleArrayGridCOrderStorageBase<T, rank, AllocationPolicy>::get(const IndexType &index) const
    {
        size_t pos = index[0];
        for (size_t i = 1; i < rank; ++i)
        {
            pos = index[i] + this->dims[i] * pos;
        }
        return this->data_fast[pos];
    }

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    inline void SingleArrayGridCOrderStorageBase<T, rank, AllocationPolicy>::resize(const IndexType &lo, const IndexType &hi) 
    {
        this->resizeImpl(lo, hi);
        // size_t p = -this->range.getLo(0);

        // for (size_t d = 1; d < rank; ++d)
        // {
        //     p = p * this->dims[d] - this->getLo(d);
        // }
        // data_fast = this->data->ptr + p;
    }

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    inline void SingleArrayGridCOrderStorageBase<T, rank, AllocationPolicy>::resize(const RangeType range)
    {
        this->resize(range.getLo(), range.getHi());
    }

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    inline ptrdiff_t SingleArrayGridCOrderStorageBase<T, rank, AllocationPolicy>::stride(size_t dim) const
    {
        size_t stride = 1;
        for (size_t i = rank - 1; i > dim; --i)
        {
            stride *= this->dims[i];
        }
        return stride;
    }

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    void SingleArrayGridCOrderStorageBase<T, rank, AllocationPolicy>::updateDataFast()
    {
        ptrdiff_t p = -this->range.getLo(0);

        for (size_t d = 1; d < rank; ++d)
        {
            p = p * this->dims[d] - this->getLo(d);
        }
        data_fast = this->data->ptr + p;
    }

    //=================================================================
    //============ SingleArrayGridFortranOrderStorageBase =============
    //=================================================================

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    SingleArrayGridFortranOrderStorageBase<T, rank, AllocationPolicy>::SingleArrayGridFortranOrderStorageBase() 
        : BaseType(), data_fast(NULL) 
    {
        this->onUpdate([this](){ updateDataFast(); });
    }

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    SingleArrayGridFortranOrderStorageBase<T, rank, AllocationPolicy>::SingleArrayGridFortranOrderStorageBase(const SingleArrayGridFortranOrderStorageBase &other)
        : BaseType(other), data_fast(other.data_fast)
    {
        this->onUpdate([this](){ updateDataFast(); });
    }

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    SingleArrayGridFortranOrderStorageBase<T, rank, AllocationPolicy>::SingleArrayGridFortranOrderStorageBase(
        const IndexType &lo, 
        const IndexType &hi
    ) : BaseType(), data_fast(NULL)
    {
        this->onUpdate([this](){ updateDataFast(); });
        resize(lo, hi);
    }

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    SingleArrayGridFortranOrderStorageBase<T, rank, AllocationPolicy>::SingleArrayGridFortranOrderStorageBase(
        const RangeType &range
    ) : BaseType(), data_fast(NULL)
    {
        this->onUpdate([this](){ updateDataFast(); });
        resize(range.getLo(), range.getHi());
    }

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    SCHNEK_INLINE T &SingleArrayGridFortranOrderStorageBase<T, rank, AllocationPolicy>::get(const IndexType &index)
    {
        size_t pos = index[rank - 1];
        for (ptrdiff_t i = ptrdiff_t(rank) - 2; i >= 0; --i)
        {
            pos = index[i] + this->dims[i] * pos;
        }
        return this->data_fast[pos];
    }

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    SCHNEK_INLINE const T &SingleArrayGridFortranOrderStorageBase<T, rank, AllocationPolicy>::get(const IndexType &index) const
    {
        size_t pos = index[rank - 1];
        for (ptrdiff_t i = ptrdiff_t(rank) - 2; i >= 0; --i)
        {
            pos = index[i] + this->dims[i] * pos;
        }
        return this->data_fast[pos];
    }

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    inline void SingleArrayGridFortranOrderStorageBase<T, rank, AllocationPolicy>::resize(const IndexType &lo, const IndexType &hi) 
    {
        this->resizeImpl(lo, hi);
    }

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    inline void SingleArrayGridFortranOrderStorageBase<T, rank, AllocationPolicy>::resize(const RangeType range)
    {
        this->resize(range.getLo(), range.getHi());
    }

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    inline ptrdiff_t SingleArrayGridFortranOrderStorageBase<T, rank, AllocationPolicy>::stride(size_t dim) const
    {
        ptrdiff_t stride = 1;
        for (size_t i = 0; i < dim; ++i)
        {
            stride *= this->dims[i];
        }
        return stride;
    }

    template <typename T, size_t rank, template <typename, size_t> class AllocationPolicy>
    void SingleArrayGridFortranOrderStorageBase<T, rank, AllocationPolicy>::updateDataFast()
    {
        size_t p = -this->getLo(rank - 1);

        for (ptrdiff_t d = ptrdiff_t(rank) - 2; d >= 0; --d)
        {
            p = p * this->dims[d] - this->getLo(d);
        }
        data_fast = this->data->ptr + p;
    }

}

#endif // SCHNEK_GRID_GRIDSTORAGE_SINGLESTORAGEBASE_HPP_