Program Listing for File MessageArray3DDevice.cuh

Return to documentation for file (include/flamegpu/runtime/messaging/MessageArray3D/MessageArray3DDevice.cuh)

#ifndef INCLUDE_FLAMEGPU_RUNTIME_MESSAGING_MESSAGEARRAY3D_MESSAGEARRAY3DDEVICE_CUH_
#define INCLUDE_FLAMEGPU_RUNTIME_MESSAGING_MESSAGEARRAY3D_MESSAGEARRAY3DDEVICE_CUH_


#include "flamegpu/runtime/messaging/MessageArray3D.h"
#include "flamegpu/runtime/messaging/MessageBruteForce/MessageBruteForceDevice.cuh"


namespace flamegpu {

class MessageArray3D::In {
    friend class Message;

 public:
    class Message {
        const MessageArray3D::In &_parent;
        size_type index;

     public:
        __device__ Message(const MessageArray3D::In &parent, const size_type _index) : _parent(parent), index(_index) {}
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
        __device__ Message(const MessageArray3D::In& parent) : _parent(parent), index(0) {}
#endif
        __device__ bool operator==(const Message& rhs) const { return  this->index == rhs.index; }
        __device__ bool operator!=(const Message& rhs) const { return  this->index != rhs.index; }
        __device__ size_type getIndex() const { return this->index; }
        template<typename T, size_type N>
        __device__ T getVariable(const char(&variable_name)[N]) const;
        template<typename T, flamegpu::size_type N, unsigned int M>
        __device__ T getVariable(const char(&variable_name)[M], unsigned int index) const;
    };
    class WrapFilter {
        friend class Message;

     public:
        class Message {
            const WrapFilter&_parent;
            int relative_cell[3];
            size_type index_1d = 0;

         public:
            __device__ Message(const WrapFilter&parent, const int relative_x, const int relative_y, const int relative_z)
                : _parent(parent) {
                relative_cell[0] = relative_x;
                relative_cell[1] = relative_y;
                relative_cell[2] = relative_z;
            }
            __device__ bool operator==(const Message& rhs) const {
                return this->relative_cell[0] == rhs.relative_cell[0]
                    && this->relative_cell[1] == rhs.relative_cell[1]
                    && this->relative_cell[2] == rhs.relative_cell[2];
                    // && this->_parent.loc[0] == rhs._parent.loc[0]
                    // && this->_parent.loc[1] == rhs._parent.loc[1]
                    // && this->_parent.loc[2] == rhs._parent.loc[2];
            }
            __device__ bool operator!=(const Message& rhs) const { return !(*this == rhs); }
            inline __device__ Message& operator++();
            __device__ size_type getX() const {
                return (this->_parent.loc[0] + relative_cell[0] + this->_parent.metadata->dimensions[0]) % this->_parent.metadata->dimensions[0];
            }
            __device__ size_type getY() const {
                return (this->_parent.loc[1] + relative_cell[1] + this->_parent.metadata->dimensions[1]) % this->_parent.metadata->dimensions[1];
            }
            __device__ size_type getZ() const {
                return (this->_parent.loc[2] + relative_cell[2] + this->_parent.metadata->dimensions[2]) % this->_parent.metadata->dimensions[2];
            }
            __device__ size_type getOffsetX() const {
                return relative_cell[0];
            }
            __device__ size_type getOffsetY() const {
                return relative_cell[1];
            }
            __device__ size_type getOffsetZ() const {
                return relative_cell[2];
            }
            template<typename T, unsigned int N>
            __device__ T getVariable(const char(&variable_name)[N]) const;
            template<typename T, flamegpu::size_type N, unsigned int M>
            __device__ T getVariable(const char(&variable_name)[M], unsigned int index) const;
        };
        class iterator {
            Message _message;

         public:
            __device__ iterator(const WrapFilter&parent, const int relative_x, const int relative_y, const int relative_z)
                : _message(parent, relative_x, relative_y, relative_z) {
                // Increment to find first message
                ++_message;
            }
            __device__ iterator& operator++() { ++_message;  return *this; }
            __device__ iterator operator++(int) {
                iterator temp = *this;
                ++*this;
                return temp;
            }
            __device__ bool operator==(const iterator& rhs) const { return  _message == rhs._message; }
            __device__ bool operator!=(const iterator& rhs) const { return  _message != rhs._message; }
            __device__ Message& operator*() { return _message; }
            __device__ Message* operator->() { return &_message; }
        };
        inline __device__ WrapFilter(const MetaData *_metadata, size_type x, size_type y, size_type z, size_type _radius);
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
        inline __device__ WrapFilter();
#endif
        inline __device__ iterator begin(void) const {
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
            if (!this->metadata)
                return iterator(*this, radius, radius, radius);
#endif
            // Bin before initial bin, as the constructor calls increment operator
            return iterator(*this, -static_cast<int>(radius), -static_cast<int>(radius), -static_cast<int>(radius)-1);
        }
        inline __device__ iterator end(void) const {
            // Final bin, as the constructor calls increment operator
            return iterator(*this, radius, radius, radius);
        }

     private:
        size_type loc[3];
        const size_type radius;
        const MetaData *metadata;
    };
    class Filter {
        friend class Message;

     public:
        class Message {
            const Filter& _parent;
            int relative_cell[3];
            size_type index_1d = 0;

         public:
            __device__ Message(const Filter& parent, const int relative_x, const int relative_y, const int relative_z)
                : _parent(parent) {
                relative_cell[0] = relative_x;
                relative_cell[1] = relative_y;
                relative_cell[2] = relative_z;
            }
            __device__ bool operator==(const Message& rhs) const {
                return this->relative_cell[0] == rhs.relative_cell[0]
                    && this->relative_cell[1] == rhs.relative_cell[1]
                    && this->relative_cell[2] == rhs.relative_cell[2];
                    // && this->_parent.loc[0] == rhs._parent.loc[0]
                    // && this->_parent.loc[1] == rhs._parent.loc[1]
                    // && this->_parent.loc[2] == rhs._parent.loc[2];
            }
            __device__ bool operator!=(const Message& rhs) const { return !(*this == rhs); }
            inline __device__ Message& operator++();
            __device__ size_type getX() const {
                return this->_parent.loc[0] + relative_cell[0];
            }
            __device__ size_type getY() const {
                return this->_parent.loc[1] + relative_cell[1];
            }
            __device__ size_type getZ() const {
                return this->_parent.loc[2] + relative_cell[2];
            }
            __device__ int getOffsetX() const {
                return this->relative_cell[0];
            }
            __device__ int getOffsetY() const {
                return this->relative_cell[1];
            }
            __device__ int getOffsetZ() const {
                return this->relative_cell[2];
            }
            template<typename T, unsigned int N>
            __device__ T getVariable(const char(&variable_name)[N]) const;
            template<typename T, flamegpu::size_type N, unsigned int M>
            __device__ T getVariable(const char(&variable_name)[M], unsigned int index) const;
        };
        class iterator {
            Message _message;

         public:
            __device__ iterator(const Filter& parent, const int relative_x, const int relative_y, const int relative_z)
                : _message(parent, relative_x, relative_y, relative_z) {
                // Increment to find first message
                ++_message;
            }
            __device__ iterator& operator++() { ++_message;  return *this; }
            __device__ iterator operator++(int) {
                iterator temp = *this;
                ++* this;
                return temp;
            }
            __device__ bool operator==(const iterator& rhs) const { return  _message == rhs._message; }
            __device__ bool operator!=(const iterator& rhs) const { return  _message != rhs._message; }
            __device__ Message& operator*() { return _message; }
            __device__ Message* operator->() { return &_message; }
        };
        inline __device__ Filter(const MetaData* _metadata, const size_type x, const size_type y, const size_type z, const size_type _radius);
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
        inline __device__ Filter();
#endif
        inline __device__ iterator begin(void) const {
            // Bin before initial bin, as the constructor calls increment operator
            return iterator(*this,  min_cell[0], min_cell[1], min_cell[2] - 1);
        }
        inline __device__ iterator end(void) const {
            // Final bin, as the constructor calls increment operator
            return iterator(*this, max_cell[0], max_cell[1], max_cell[2]);
        }

     private:
        size_type loc[3];
        int min_cell[3];
        int max_cell[3];
        const MetaData* metadata;
    };
    class VonNeumannWrapFilter {
        friend class Message;

     public:
        class Message {
            const VonNeumannWrapFilter&_parent;
            int relative_cell[3];
            size_type index_1d = 0;

         public:
            __device__ Message(const VonNeumannWrapFilter&parent, const int relative_x, const int relative_y, const int relative_z)
                : _parent(parent) {
                relative_cell[0] = relative_x;
                relative_cell[1] = relative_y;
                relative_cell[2] = relative_z;
            }
            __device__ bool operator==(const Message& rhs) const {
                return this->relative_cell[0] == rhs.relative_cell[0]
                    && this->relative_cell[1] == rhs.relative_cell[1]
                    && this->relative_cell[2] == rhs.relative_cell[2];
                    // && this->_parent.loc[0] == rhs._parent.loc[0]
                    // && this->_parent.loc[1] == rhs._parent.loc[1]
                    // && this->_parent.loc[2] == rhs._parent.loc[2];
            }
            __device__ bool operator!=(const Message& rhs) const { return !(*this == rhs); }
            inline __device__ Message& operator++();
            __device__ size_type getX() const {
                return (this->_parent.loc[0] + relative_cell[0] + this->_parent.metadata->dimensions[0]) % this->_parent.metadata->dimensions[0];
            }
            __device__ size_type getY() const {
                return (this->_parent.loc[1] + relative_cell[1] + this->_parent.metadata->dimensions[1]) % this->_parent.metadata->dimensions[1];
            }
            __device__ size_type getZ() const {
                return (this->_parent.loc[2] + relative_cell[2] + this->_parent.metadata->dimensions[2]) % this->_parent.metadata->dimensions[2];
            }
            __device__ size_type getOffsetX() const {
                return relative_cell[0];
            }
            __device__ size_type getOffsetY() const {
                return relative_cell[1];
            }
            __device__ size_type getOffsetZ() const {
                return relative_cell[2];
            }
            template<typename T, unsigned int N>
            __device__ T getVariable(const char(&variable_name)[N]) const;
            template<typename T, flamegpu::size_type N, unsigned int M>
            __device__ T getVariable(const char(&variable_name)[M], unsigned int index) const;
        };
        class iterator {
            Message _message;

         public:
            __device__ iterator(const VonNeumannWrapFilter&parent, const int relative_x, const int relative_y, const int relative_z)
                : _message(parent, relative_x, relative_y, relative_z) {
                // Increment to find first message
                ++_message;
            }
            __device__ iterator& operator++() { ++_message;  return *this; }
            __device__ iterator operator++(int) {
                iterator temp = *this;
                ++*this;
                return temp;
            }
            __device__ bool operator==(const iterator& rhs) const { return  _message == rhs._message; }
            __device__ bool operator!=(const iterator& rhs) const { return  _message != rhs._message; }
            __device__ Message& operator*() { return _message; }
            __device__ Message* operator->() { return &_message; }
        };
        inline __device__ VonNeumannWrapFilter(const MetaData *_metadata, size_type x, size_type y, size_type z, size_type _radius);
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
        inline __device__ VonNeumannWrapFilter();
#endif
        inline __device__ iterator begin(void) const {
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
            if (!this->metadata)
                return iterator(*this, radius, radius, radius);
#endif
            // Bin before initial bin, as the constructor calls increment operator
            return iterator(*this, -static_cast<int>(radius), -static_cast<int>(radius), -static_cast<int>(radius)-1);
        }
        inline __device__ iterator end(void) const {
            // Final bin, as the constructor calls increment operator
            return iterator(*this, radius, radius, radius);
        }

     private:
        size_type loc[3];
        const int radius;
        const MetaData *metadata;
    };
    class VonNeumannFilter {
        friend class Message;

     public:
        class Message {
            const VonNeumannFilter& _parent;
            int relative_cell[3];
            size_type index_1d = 0;

         public:
            __device__ Message(const VonNeumannFilter& parent, const int relative_x, const int relative_y, const int relative_z)
                : _parent(parent) {
                relative_cell[0] = relative_x;
                relative_cell[1] = relative_y;
                relative_cell[2] = relative_z;
            }
            __device__ bool operator==(const Message& rhs) const {
                return this->relative_cell[0] == rhs.relative_cell[0]
                    && this->relative_cell[1] == rhs.relative_cell[1]
                    && this->relative_cell[2] == rhs.relative_cell[2];
                    // && this->_parent.loc[0] == rhs._parent.loc[0]
                    // && this->_parent.loc[1] == rhs._parent.loc[1]
                    // && this->_parent.loc[2] == rhs._parent.loc[2];
            }
            __device__ bool operator!=(const Message& rhs) const { return !(*this == rhs); }
            inline __device__ Message& operator++();
            __device__ size_type getX() const {
                return this->_parent.loc[0] + relative_cell[0];
            }
            __device__ size_type getY() const {
                return this->_parent.loc[1] + relative_cell[1];
            }
            __device__ size_type getZ() const {
                return this->_parent.loc[2] + relative_cell[2];
            }
            __device__ int getOffsetX() const {
                return this->relative_cell[0];
            }
            __device__ int getOffsetY() const {
                return this->relative_cell[1];
            }
            __device__ int getOffsetZ() const {
                return this->relative_cell[2];
            }
            template<typename T, unsigned int N>
            __device__ T getVariable(const char(&variable_name)[N]) const;
            template<typename T, flamegpu::size_type N, unsigned int M>
            __device__ T getVariable(const char(&variable_name)[M], unsigned int index) const;
        };
        class iterator {
            Message _message;

         public:
            __device__ iterator(const VonNeumannFilter& parent, const int relative_x, const int relative_y, const int relative_z)
                : _message(parent, relative_x, relative_y, relative_z) {
                // Increment to find first message
                ++_message;
            }
            __device__ iterator& operator++() { ++_message;  return *this; }
            __device__ iterator operator++(int) {
                iterator temp = *this;
                ++* this;
                return temp;
            }
            __device__ bool operator==(const iterator& rhs) const { return  _message == rhs._message; }
            __device__ bool operator!=(const iterator& rhs) const { return  _message != rhs._message; }
            __device__ Message& operator*() { return _message; }
            __device__ Message* operator->() { return &_message; }
        };
        inline __device__ VonNeumannFilter(const MetaData* _metadata, size_type x, size_type y, size_type z, size_type _radius);
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
        inline __device__ VonNeumannFilter();
#endif
        inline __device__ iterator begin(void) const {
            // Bin before initial bin, as the constructor calls increment operator
            return iterator(*this,  min_cell[0], min_cell[1], min_cell[2] - 1);
        }
        inline __device__ iterator end(void) const {
            // Final bin, as the constructor calls increment operator
            return iterator(*this, max_cell[0], max_cell[1], max_cell[2]);
        }

     private:
        size_type loc[3];
        int min_cell[3];
        int max_cell[3];
        const int radius;
        const MetaData* metadata;
    };
    __device__ In(const void *_metadata)
        : metadata(reinterpret_cast<const MetaData*>(_metadata))
    { }
    inline __device__ WrapFilter wrap(const size_type x, const size_type y, const size_type z, const size_type radius = 1) const {
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
        if (radius == 0) {
            DTHROW("%u is not a valid radius for accessing Array3D message lists.\n", radius);
            return WrapFilter();
        } else if ((radius * 2) + 1 > metadata->dimensions[0] ||
                   (radius * 2) + 1 > metadata->dimensions[1] ||
                   (radius * 2) + 1 > metadata->dimensions[2]) {
            unsigned int min_r = metadata->dimensions[0] < metadata->dimensions[1] ? metadata->dimensions[0] : metadata->dimensions[1];
            min_r = min_r < metadata->dimensions[2] ? min_r : metadata->dimensions[2];
            min_r = min_r % 2 == 0 ? min_r - 2: min_r - 1;
            min_r /= 2;
            DTHROW("%u is not a valid radius for accessing Array3D message lists, as the diameter of messages accessed exceeds one or more of the message list dimensions (%u, %u, %u)."
            " Maximum supported radius for this message list is %u.\n",
            radius, metadata->dimensions[0], metadata->dimensions[1], metadata->dimensions[2], min_r);
            return WrapFilter();
        } else if (x >= metadata->dimensions[0] ||
                   y >= metadata->dimensions[1] ||
                   z >= metadata->dimensions[2]) {
            DTHROW("(%u, %u, %u) is not a valid position for iterating an Array3D message list of dimensions (%u, %u, %u), location must be within bounds.",
                x, y, z, metadata->dimensions[0], metadata->dimensions[1], metadata->dimensions[2]);
            return WrapFilter();
        }
#endif
        return WrapFilter(metadata, x, y, z, radius);
    }
    inline __device__ Filter operator()(const size_type x, const size_type y, const size_type z, const size_type radius = 1) const {
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
        if (radius == 0) {
            DTHROW("%u is not a valid radius for accessing Array3D message lists.\n", radius);
            return Filter();
        } else if (x >= metadata->dimensions[0] ||
                   y >= metadata->dimensions[1] ||
                   z >= metadata->dimensions[2]) {
            DTHROW("(%u, %u, %u) is not a valid position for iterating an Array3D message list of dimensions (%u, %u, %u), location must be within bounds.",
                x, y, z, metadata->dimensions[0], metadata->dimensions[1], metadata->dimensions[2]);
            return Filter();
        }
#endif
        return Filter(metadata, x, y, z, radius);
    }
    inline __device__ VonNeumannWrapFilter vn_wrap(const size_type x, const size_type y, const size_type z, const size_type radius = 1) const {
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
        if (radius == 0) {
            DTHROW("%u is not a valid radius for accessing Array3D message lists.\n", radius);
            return VonNeumannWrapFilter();
        } else if ((radius * 2) + 1 > metadata->dimensions[0] ||
                   (radius * 2) + 1 > metadata->dimensions[1] ||
                   (radius * 2) + 1 > metadata->dimensions[2]) {
            unsigned int min_r = metadata->dimensions[0] < metadata->dimensions[1] ? metadata->dimensions[0] : metadata->dimensions[1];
            min_r = min_r < metadata->dimensions[2] ? min_r : metadata->dimensions[2];
            min_r = min_r % 2 == 0 ? min_r - 2: min_r - 1;
            min_r /= 2;
            DTHROW("%u is not a valid radius for accessing Array3D message lists, as the diameter of messages accessed exceeds one or more of the message list dimensions (%u, %u, %u)."
            " Maximum supported radius for this message list is %u.\n",
            radius, metadata->dimensions[0], metadata->dimensions[1], metadata->dimensions[2], min_r);
            return VonNeumannWrapFilter();
        } else if (x >= metadata->dimensions[0] ||
                   y >= metadata->dimensions[1] ||
                   z >= metadata->dimensions[2]) {
            DTHROW("(%u, %u, %u) is not a valid position for iterating an Array3D message list of dimensions (%u, %u, %u), location must be within bounds.",
                x, y, z, metadata->dimensions[0], metadata->dimensions[1], metadata->dimensions[2]);
            return VonNeumannWrapFilter();
        }
#endif
        return VonNeumannWrapFilter(metadata, x, y, z, radius);
    }
    inline __device__ VonNeumannFilter vn(const size_type x, const size_type y, const size_type z, const size_type radius = 1) const {
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
        if (radius == 0) {
            DTHROW("%u is not a valid radius for accessing Array3D message lists.\n", radius);
            return VonNeumannFilter();
        } else if (x >= metadata->dimensions[0] ||
                   y >= metadata->dimensions[1] ||
                   z >= metadata->dimensions[2]) {
            DTHROW("(%u, %u, %u) is not a valid position for iterating an Array3D message list of dimensions (%u, %u, %u), location must be within bounds.",
                x, y, z, metadata->dimensions[0], metadata->dimensions[1], metadata->dimensions[2]);
            return VonNeumannFilter();
        }
#endif
        return VonNeumannFilter(metadata, x, y, z, radius);
    }
    __device__ size_type getDimX() const {
        return metadata->dimensions[0];
    }
    __device__ size_type getDimY() const {
        return metadata->dimensions[1];
    }
    __device__ size_type getDimZ() const {
        return metadata->dimensions[2];
    }
    __device__ size_type size(void) const {
        return metadata->length;
    }
    __device__ Message at(const size_type x, const size_type y, const size_type z) const {
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
        if (x >= metadata->dimensions[0] || y >= metadata->dimensions[1] || z >= metadata->dimensions[2]) {
            DTHROW("Index is out of bounds for Array3D messagelist ([%u, %u, %u] >= [%u, %u, %u]).\n", x, y, z, metadata->dimensions[0], metadata->dimensions[1], metadata->dimensions[2]);
            return Message(*this);
        }
#endif
        const size_type index_1d =
            z * metadata->dimensions[0] * metadata->dimensions[1]  +
            y * metadata->dimensions[0]  +
            x;
        return Message(*this, index_1d);
    }

 private:
    const MetaData * const metadata;
};

class MessageArray3D::Out {
 public:
    __device__ Out(const void *_metadata, unsigned int *scan_flag_messageOutput)
        : scan_flag(scan_flag_messageOutput)
        , metadata(reinterpret_cast<const MetaData*>(_metadata))
    { }
    inline __device__ void setIndex(size_type x, size_type y, size_type z) const;
    template<typename T, unsigned int N>
    __device__ void setVariable(const char(&variable_name)[N], T value) const;
    template<typename T, unsigned int N, unsigned int M>
    __device__ void setVariable(const char(&variable_name)[M], unsigned int index, T value) const;

 protected:
    unsigned int *scan_flag;
    const MetaData * const metadata;
};

template<typename T, unsigned int N>
__device__ T MessageArray3D::In::Message::getVariable(const char(&variable_name)[N]) const {
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
    // Ensure that the message is within bounds.
    if (index >= this->_parent.metadata->length) {
        DTHROW("Invalid Array3D message, unable to get variable '%s'.\n", variable_name);
        return static_cast<T>(0);
    }
#endif
    // get the value from curve using the message index.
    return detail::curve::DeviceCurve::getMessageVariable<T>(variable_name, index);
}
template<typename T, flamegpu::size_type N, unsigned int M> __device__
T MessageArray3D::In::Message::getVariable(const char(&variable_name)[M], const unsigned int array_index) const {
    // simple indexing assumes index is the thread number (this may change later)
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
    // Ensure that the message is within bounds.
    if (index >= this->_parent.metadata->length) {
        DTHROW("Invalid Array3D message, unable to get variable '%s'.\n", variable_name);
        return {};
    }
#endif
    // get the value from curve using the message index.
    T value = detail::curve::DeviceCurve::getMessageArrayVariable<T, N>(variable_name, index, array_index);
    return value;
}
template<typename T, unsigned int N>
__device__ T MessageArray3D::In::WrapFilter::Message::getVariable(const char(&variable_name)[N]) const {
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
    // Ensure that the message is within bounds.
    if (index_1d >= this->_parent.metadata->length) {
        DTHROW("Invalid Array3D message, unable to get variable '%s'.\n", variable_name);
        return static_cast<T>(0);
    }
#endif
    // get the value from curve using the message index.
    return detail::curve::DeviceCurve::getMessageVariable<T>(variable_name, index_1d);
}
template<typename T, flamegpu::size_type N, unsigned int M> __device__
T MessageArray3D::In::WrapFilter::Message::getVariable(const char(&variable_name)[M], const unsigned int array_index) const {
    // simple indexing assumes index is the thread number (this may change later)
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
    // Ensure that the message is within bounds.
    if (index_1d >= this->_parent.metadata->length) {
        DTHROW("Invalid Array3D message, unable to get variable '%s'.\n", variable_name);
        return {};
    }
#endif
    // get the value from curve using the message index.
    T value = detail::curve::DeviceCurve::getMessageArrayVariable<T, N>(variable_name, index_1d, array_index);
    return value;
}
template<typename T, unsigned int N>
__device__ T MessageArray3D::In::Filter::Message::getVariable(const char(&variable_name)[N]) const {
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
    // Ensure that the message is within bounds.
    if (index_1d >= this->_parent.metadata->length) {
        DTHROW("Invalid Array3D message, unable to get variable '%s'.\n", variable_name);
        return static_cast<T>(0);
    }
#endif
    // get the value from curve using the message index.
    return detail::curve::DeviceCurve::getMessageVariable<T>(variable_name, index_1d);
}
template<typename T, flamegpu::size_type N, unsigned int M> __device__
T MessageArray3D::In::Filter::Message::getVariable(const char(&variable_name)[M], const unsigned int array_index) const {
    // simple indexing assumes index is the thread number (this may change later)
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
    // Ensure that the message is within bounds.
    if (index_1d >= this->_parent.metadata->length) {
        DTHROW("Invalid Array3D message, unable to get variable '%s'.\n", variable_name);
        return {};
    }
#endif
    // get the value from curve using the message index.
    T value = detail::curve::DeviceCurve::getMessageArrayVariable<T, N>(variable_name, index_1d, array_index);
    return value;
}
template<typename T, unsigned int N>
__device__ T MessageArray3D::In::VonNeumannWrapFilter::Message::getVariable(const char(&variable_name)[N]) const {
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
    // Ensure that the message is within bounds.
    if (index_1d >= this->_parent.metadata->length) {
        DTHROW("Invalid Array3D message, unable to get variable '%s'.\n", variable_name);
        return static_cast<T>(0);
    }
#endif
    // get the value from curve using the message index.
    return detail::curve::DeviceCurve::getMessageVariable<T>(variable_name, index_1d);
}
template<typename T, flamegpu::size_type N, unsigned int M> __device__
T MessageArray3D::In::VonNeumannWrapFilter::Message::getVariable(const char(&variable_name)[M], const unsigned int array_index) const {
    // simple indexing assumes index is the thread number (this may change later)
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
    // Ensure that the message is within bounds.
    if (index_1d >= this->_parent.metadata->length) {
        DTHROW("Invalid Array3D message, unable to get variable '%s'.\n", variable_name);
        return {};
    }
#endif
    // get the value from curve using the message index.
    T value = detail::curve::DeviceCurve::getMessageArrayVariable<T, N>(variable_name, index_1d, array_index);
    return value;
}
template<typename T, unsigned int N>
__device__ T MessageArray3D::In::VonNeumannFilter::Message::getVariable(const char(&variable_name)[N]) const {
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
    // Ensure that the message is within bounds.
    if (index_1d >= this->_parent.metadata->length) {
        DTHROW("Invalid Array3D message, unable to get variable '%s'.\n", variable_name);
        return static_cast<T>(0);
    }
#endif
    // get the value from curve using the message index.
    return detail::curve::DeviceCurve::getMessageVariable<T>(variable_name, index_1d);
}
template<typename T, flamegpu::size_type N, unsigned int M> __device__
T MessageArray3D::In::VonNeumannFilter::Message::getVariable(const char(&variable_name)[M], const unsigned int array_index) const {
    // simple indexing assumes index is the thread number (this may change later)
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
    // Ensure that the message is within bounds.
    if (index_1d >= this->_parent.metadata->length) {
        DTHROW("Invalid Array3D message, unable to get variable '%s'.\n", variable_name);
        return {};
    }
#endif
    // get the value from curve using the message index.
    T value = detail::curve::DeviceCurve::getMessageArrayVariable<T, N>(variable_name, index_1d, array_index);
    return value;
}

template<typename T, unsigned int N>
__device__ void MessageArray3D::Out::setVariable(const char(&variable_name)[N], T value) const {  // message name or variable name
    if (variable_name[0] == '_') {
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
        DTHROW("Variable names starting with '_' are reserved for internal use, with '%s', in MessageArray3D::Out::setVariable().\n", variable_name);
#endif
        return;  // Fail silently
    }
    unsigned int index = (blockDim.x * blockIdx.x) + threadIdx.x;

    // set the variable using curve
    detail::curve::DeviceCurve::setMessageVariable<T>(variable_name, value, index);

    // setIndex() sets the optional message scan flag
}
template<typename T, unsigned int N, unsigned int M>
__device__ void MessageArray3D::Out::setVariable(const char(&variable_name)[M], const unsigned int array_index, T value) const {
    if (variable_name[0] == '_') {
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
        DTHROW("Variable names starting with '_' are reserved for internal use, with '%s', in MessageArray3D::Out::setVariable().\n", variable_name);
#endif
        return;  // Fail silently
    }
    unsigned int index = (blockDim.x * blockIdx.x) + threadIdx.x;

    // set the variable using curve
    detail::curve::DeviceCurve::setMessageArrayVariable<T, N>(variable_name, value, index, array_index);

    // setIndex() sets the optional message scan flag
}

__device__ inline void MessageArray3D::Out::setIndex(const size_type x, const size_type y, const size_type z) const {
    unsigned int index = (blockDim.x * blockIdx.x) + threadIdx.x;
    size_type index_1d =
        z * metadata->dimensions[0] * metadata->dimensions[1] +
        y * metadata->dimensions[0] +
        x;
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
    if (x >= metadata->dimensions[0] ||
        y >= metadata->dimensions[1] ||
        z >= metadata->dimensions[2]) {
        DTHROW("MessageArray3D index [%u, %u, %u] is out of bounds [%u, %u, %u]\n", x, y, z, metadata->dimensions[0], metadata->dimensions[1], metadata->dimensions[2]);
        return;  // Fail silently
    }
#endif

    // set the variable using curve
    detail::curve::DeviceCurve::setMessageVariable<size_type>("___INDEX", index_1d, index);

    // Set scan flag incase the message is optional
    this->scan_flag[index] = 1;
}
// Moore Wrap
__device__ inline MessageArray3D::In::WrapFilter::WrapFilter(const MetaData *_metadata, const size_type x, const size_type y, const size_type z, const size_type _radius)
    : radius(_radius)
    , metadata(_metadata) {
    loc[0] = x;
    loc[1] = y;
    loc[2] = z;
}
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
__device__ inline MessageArray3D::In::WrapFilter::WrapFilter()
    : radius(0)
    , metadata(nullptr) {
    loc[0] = 0;
    loc[1] = 0;
    loc[2] = 0;
}
#endif
__device__ inline MessageArray3D::In::WrapFilter::Message& MessageArray3D::In::WrapFilter::Message::operator++() {
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
    if (!_parent.metadata)
        return *this;
#endif
    if (relative_cell[2] >= static_cast<int>(_parent.radius)) {
        relative_cell[2] = -static_cast<int>(_parent.radius);
        if (relative_cell[1] >= static_cast<int>(_parent.radius)) {
            relative_cell[1] = -static_cast<int>(_parent.radius);
            relative_cell[0]++;
        } else {
            relative_cell[1]++;
        }
    } else {
        relative_cell[2]++;
    }
    // Skip origin cell
    if (relative_cell[0] == 0 && relative_cell[1] == 0 && relative_cell[2] == 0) {
        relative_cell[2]++;
    }
    // Wrap over boundaries
    const unsigned int their_x = (this->_parent.loc[0] + relative_cell[0] + this->_parent.metadata->dimensions[0]) % this->_parent.metadata->dimensions[0];
    const unsigned int their_y = (this->_parent.loc[1] + relative_cell[1] + this->_parent.metadata->dimensions[1]) % this->_parent.metadata->dimensions[1];
    const unsigned int their_z = (this->_parent.loc[2] + relative_cell[2] + this->_parent.metadata->dimensions[2]) % this->_parent.metadata->dimensions[2];
    // Solve to 1 dimensional bin index
    index_1d = their_z * this->_parent.metadata->dimensions[0] * this->_parent.metadata->dimensions[1] +
               their_y * this->_parent.metadata->dimensions[0] +
               their_x;
    return *this;
}
// Moore
__device__ inline MessageArray3D::In::Filter::Filter(const MetaData* _metadata, const size_type x, const size_type y, const size_type z, const size_type _radius)
    : metadata(_metadata) {
    loc[0] = x;
    loc[1] = y;
    loc[2] = z;
    min_cell[0] = static_cast<int>(x) - static_cast<int>(_radius) < 0 ? -static_cast<int>(x) : - static_cast<int>(_radius);
    min_cell[1] = static_cast<int>(y) - static_cast<int>(_radius) < 0 ? -static_cast<int>(y) : - static_cast<int>(_radius);
    min_cell[2] = static_cast<int>(z) - static_cast<int>(_radius) < 0 ? -static_cast<int>(z) : - static_cast<int>(_radius);
    max_cell[0] = x + _radius >= _metadata->dimensions[0] ? static_cast<int>(_metadata->dimensions[0]) - 1 - static_cast<int>(x) : static_cast<int>(_radius);
    max_cell[1] = y + _radius >= _metadata->dimensions[1] ? static_cast<int>(_metadata->dimensions[1]) - 1 - static_cast<int>(y) : static_cast<int>(_radius);
    max_cell[2] = z + _radius >= _metadata->dimensions[2] ? static_cast<int>(_metadata->dimensions[2]) - 1 - static_cast<int>(z) : static_cast<int>(_radius);
}
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
__device__ inline MessageArray3D::In::Filter::Filter()
    : metadata(nullptr) {
    loc[0] = 0;
    loc[1] = 0;
    loc[2] = 0;
    min_cell[0] = 0;
    min_cell[1] = 0;
    min_cell[2] = 1;
    max_cell[0] = 0;
    max_cell[1] = 0;
    max_cell[2] = 0;
}
#endif
__device__ inline MessageArray3D::In::Filter::Message& MessageArray3D::In::Filter::Message::operator++() {
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
    if (!_parent.metadata)
        return *this;
#endif
    if (relative_cell[2] >= _parent.max_cell[2]) {
        relative_cell[2] = _parent.min_cell[2];
        if (relative_cell[1] >= _parent.max_cell[1]) {
            relative_cell[1] = _parent.min_cell[1];
            relative_cell[0]++;
        } else {
            relative_cell[1]++;
        }
    } else {
        relative_cell[2]++;
    }
    // Skip origin cell
    if (relative_cell[0] == 0 && relative_cell[1] == 0 && relative_cell[2] == 0) {
        if (relative_cell[2] >= _parent.max_cell[2]) {
            relative_cell[2] = _parent.min_cell[2];
            if (relative_cell[1] >= _parent.max_cell[1]) {
                relative_cell[1] = _parent.min_cell[1];
                relative_cell[0]++;
            } else {
                relative_cell[1]++;
            }
        } else {
            relative_cell[2]++;
        }
    }
    // Solve to 1 dimensional bin index
    index_1d = (this->_parent.loc[2] + relative_cell[2]) * this->_parent.metadata->dimensions[0] * this->_parent.metadata->dimensions[1] +
               (this->_parent.loc[1] + relative_cell[1]) * this->_parent.metadata->dimensions[0] +
               (this->_parent.loc[0] + relative_cell[0]);
    return *this;
}
// Von Neumann Wrap
__device__ inline MessageArray3D::In::VonNeumannWrapFilter::VonNeumannWrapFilter(const MetaData *_metadata, const size_type x, const size_type y, const size_type z, const size_type _radius)
    : radius(static_cast<int>(_radius))
    , metadata(_metadata) {
    loc[0] = x;
    loc[1] = y;
    loc[2] = z;
}
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
__device__ inline MessageArray3D::In::VonNeumannWrapFilter::VonNeumannWrapFilter()
    : radius(0)
    , metadata(nullptr) {
    loc[0] = 0;
    loc[1] = 0;
    loc[2] = 0;
}
#endif
__device__ inline MessageArray3D::In::VonNeumannWrapFilter::Message& MessageArray3D::In::VonNeumannWrapFilter::Message::operator++() {
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
    if (!_parent.metadata)
        return *this;
#endif
    do {
        if (relative_cell[2] >= static_cast<int>(_parent.radius)) {
            relative_cell[2] = -static_cast<int>(_parent.radius);
            if (relative_cell[1] >= static_cast<int>(_parent.radius)) {
                relative_cell[1] = -static_cast<int>(_parent.radius);
                relative_cell[0]++;
            } else {
                relative_cell[1]++;
            }
        } else {
            relative_cell[2]++;
        }
        // Skip origin cell
        if (relative_cell[0] == 0 && relative_cell[1] == 0 && relative_cell[2] == 0) {
            relative_cell[2]++;
        }
        // Break if we reach the end of the Moore neighbourhood
        if (relative_cell[0] >= static_cast<int>(_parent.radius) + 1) {
            break;
        }
    // Skip cells outside the valid manhattan distance
    } while (abs(relative_cell[0]) + abs(relative_cell[1]) + abs(relative_cell[2]) > static_cast<int>(_parent.radius));
    // Wrap over boundaries
    const unsigned int their_x = (this->_parent.loc[0] + relative_cell[0] + this->_parent.metadata->dimensions[0]) % this->_parent.metadata->dimensions[0];
    const unsigned int their_y = (this->_parent.loc[1] + relative_cell[1] + this->_parent.metadata->dimensions[1]) % this->_parent.metadata->dimensions[1];
    const unsigned int their_z = (this->_parent.loc[2] + relative_cell[2] + this->_parent.metadata->dimensions[2]) % this->_parent.metadata->dimensions[2];
    // Solve to 1 dimensional bin index
    index_1d = their_z * this->_parent.metadata->dimensions[0] * this->_parent.metadata->dimensions[1] +
               their_y * this->_parent.metadata->dimensions[0] +
               their_x;
    return *this;
}
// Von Neumann
__device__ inline MessageArray3D::In::VonNeumannFilter::VonNeumannFilter(const MetaData* _metadata, const size_type x, const size_type y, const size_type z, const size_type _radius)
    : radius(static_cast<int>(_radius))
    , metadata(_metadata) {
    loc[0] = x;
    loc[1] = y;
    loc[2] = z;
    min_cell[0] = static_cast<int>(x) - static_cast<int>(_radius) < 0 ? -static_cast<int>(x) : - static_cast<int>(_radius);
    min_cell[1] = static_cast<int>(y) - static_cast<int>(_radius) < 0 ? -static_cast<int>(y) : - static_cast<int>(_radius);
    min_cell[2] = static_cast<int>(z) - static_cast<int>(_radius) < 0 ? -static_cast<int>(z) : - static_cast<int>(_radius);
    max_cell[0] = x + _radius >= _metadata->dimensions[0] ? static_cast<int>(_metadata->dimensions[0]) - 1 - static_cast<int>(x) : static_cast<int>(_radius);
    max_cell[1] = y + _radius >= _metadata->dimensions[1] ? static_cast<int>(_metadata->dimensions[1]) - 1 - static_cast<int>(y) : static_cast<int>(_radius);
    max_cell[2] = z + _radius >= _metadata->dimensions[2] ? static_cast<int>(_metadata->dimensions[2]) - 1 - static_cast<int>(z) : static_cast<int>(_radius);
}
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
__device__ inline MessageArray3D::In::VonNeumannFilter::VonNeumannFilter()
    : radius(0)
    , metadata(nullptr) {
    loc[0] = 0;
    loc[1] = 0;
    loc[2] = 0;
    min_cell[0] = 0;
    min_cell[1] = 0;
    min_cell[2] = 1;
    max_cell[0] = 0;
    max_cell[1] = 0;
    max_cell[2] = 0;
}
#endif
__device__ inline MessageArray3D::In::VonNeumannFilter::Message& MessageArray3D::In::VonNeumannFilter::Message::operator++() {
#if !defined(FLAMEGPU_SEATBELTS) || FLAMEGPU_SEATBELTS
    if (!_parent.metadata)
        return *this;
#endif
    do {
        if (relative_cell[2] >= _parent.max_cell[2]) {
            relative_cell[2] = _parent.min_cell[2];
            if (relative_cell[1] >= _parent.max_cell[1]) {
                relative_cell[1] = _parent.min_cell[1];
                relative_cell[0]++;
            } else {
                relative_cell[1]++;
            }
        } else {
            relative_cell[2]++;
        }
        // Skip origin cell
        if (relative_cell[0] == 0 && relative_cell[1] == 0 && relative_cell[2] == 0) {
            if (relative_cell[2] >= _parent.max_cell[2]) {
                relative_cell[2] = _parent.min_cell[2];
                if (relative_cell[1] >= _parent.max_cell[1]) {
                    relative_cell[1] = _parent.min_cell[1];
                    relative_cell[0]++;
                } else {
                    relative_cell[1]++;
                }
            } else {
                relative_cell[2]++;
            }
        }
        // Break if we reach the end of the Moore neighbourhood
        if (relative_cell[0] >= _parent.max_cell[0] + 1) {
            break;
        }
    // Skip cells outside the valid manhattan distance
    } while (abs(relative_cell[0]) + abs(relative_cell[1]) + abs(relative_cell[2]) > static_cast<int>(_parent.radius));
    // Solve to 1 dimensional bin index
    index_1d = (this->_parent.loc[2] + relative_cell[2]) * this->_parent.metadata->dimensions[0] * this->_parent.metadata->dimensions[1] +
               (this->_parent.loc[1] + relative_cell[1]) * this->_parent.metadata->dimensions[0] +
               (this->_parent.loc[0] + relative_cell[0]);
    return *this;
}

}  // namespace flamegpu

#endif  // INCLUDE_FLAMEGPU_RUNTIME_MESSAGING_MESSAGEARRAY3D_MESSAGEARRAY3DDEVICE_CUH_