cube.hpp

Go to the documentation of this file.

/*
Copyright (c) 2014-2016, Scott Zuyderduyn
All rights reserved.

Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions are met:

1. Redistributions of source code must retain the above copyright notice, this
   list of conditions and the following disclaimer.
2. Redistributions in binary form must reproduce the above copyright notice,
   this list of conditions and the following disclaimer in the documentation
   and/or other materials provided with the distribution.

THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR
ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
(INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

The views and conclusions contained in the software and documentation are those
of the authors and should not be interpreted as representing official policies,
either expressed or implied, of the FreeBSD Project.
*/

//----------------------------------------------------------------------------
// cube.hpp
// An STL-like structure that resides in video memory.
//
// Author: Scott D. Zuyderduyn, Ph.D. (scott.zuyderduyn@utoronto.ca)
//----------------------------------------------------------------------------

#pragma once
#ifndef ECUDA_CUBE_HPP
#define ECUDA_CUBE_HPP

#include <limits>
#include <stdexcept>
#include <vector>

#include "global.hpp"
#include "allocators.hpp"
#include "matrix.hpp"
#include "memory.hpp"
#include "type_traits.hpp"

#include "model/device_sequence.hpp"
#include "model/device_contiguous_sequence.hpp"
#include "model/device_matrix.hpp"
#include "model/device_contiguous_row_matrix.hpp"

namespace ecuda {

namespace impl {

template<typename T,class Alloc> class cube_kernel_argument; // forward declaration

} // namespace impl

template< typename T, class Alloc=device_pitch_allocator<T>, class P=shared_ptr<T> >
class cube : private model::device_contiguous_row_matrix< T, /*padded_ptr< T,*/ P /*>*/ > {

private:
    typedef model::device_contiguous_row_matrix< T, /*padded_ptr< T,*/ P /*>*/ > base_type;

public:
    typedef typename base_type::value_type      value_type;      
    typedef Alloc                               allocator_type;  
    typedef typename base_type::size_type       size_type;       
    typedef typename base_type::difference_type difference_type; 
    typedef typename base_type::reference       reference;       
    typedef typename base_type::const_reference const_reference; 
    typedef typename base_type::pointer         pointer;         
    typedef typename make_const<pointer>::type  const_pointer;   

    typedef model::device_sequence<           value_type, striding_padded_ptr<value_type,typename ecuda::add_pointer<value_type>::type> > row_type;    
    typedef model::device_sequence<           value_type, striding_padded_ptr<value_type,typename ecuda::add_pointer<value_type>::type> > column_type; 
    typedef model::device_contiguous_sequence<value_type                                                                                > depth_type;  
    typedef model::device_sequence<           const value_type, striding_padded_ptr<const value_type,typename ecuda::add_pointer<const value_type>::type> > const_row_type;    
    typedef model::device_sequence<           const value_type, striding_padded_ptr<const value_type,typename ecuda::add_pointer<const value_type>::type> > const_column_type; 
    typedef model::device_contiguous_sequence<const value_type                                                                                            > const_depth_type;  

    typedef typename base_type::iterator               iterator;               
    typedef typename base_type::const_iterator         const_iterator;         
    typedef typename base_type::reverse_iterator       reverse_iterator;       
    typedef typename base_type::const_reverse_iterator const_reverse_iterator; 

    typedef model::device_matrix<                value_type,       striding_padded_ptr<value_type,typename ecuda::add_pointer<value_type>::type> > slice_xy_type; 
    typedef model::device_contiguous_row_matrix< value_type,       typename ecuda::add_pointer<value_type>::type                                 > slice_xz_type; 
    typedef model::device_contiguous_row_matrix< value_type,       typename ecuda::add_pointer<value_type>::type                                 > slice_yz_type; 
    typedef model::device_matrix<                const value_type, striding_padded_ptr<const value_type,typename ecuda::add_pointer<const value_type>::type> > const_slice_xy_type; 
    typedef model::device_contiguous_row_matrix< const value_type, typename ecuda::add_pointer<const value_type>::type                                       > const_slice_xz_type; 
    typedef model::device_contiguous_row_matrix< const value_type, typename ecuda::add_pointer<const value_type>::type                                       > const_slice_yz_type; 

    typedef       impl::cube_kernel_argument<T,Alloc> kernel_argument;       
    typedef const impl::cube_kernel_argument<T,Alloc> const_kernel_argument; 

private:
    size_type numberRows; 
    allocator_type allocator; 

    template<typename U,class Alloc2,class Q> friend class cube;

protected:
    template<typename U>
    __HOST__ __DEVICE__ cube( const cube<T,Alloc,U>& src, ecuda::true_type ) : base_type( unmanaged_cast(src.get_pointer()), src.number_rows()*src.number_columns(), src.number_depths() ), numberRows(src.numberRows), allocator(src.allocator) {}

    __HOST__ __DEVICE__ cube& shallow_assign( const cube& other )
    {
        base_type::get_pointer() = other.get_pointer();
        allocator = other.allocator;
        numberRows = other.numberRows;
        return *this;
    }

private:
    __HOST__ void init()
    {
        if( number_rows() && number_columns() && number_depths() ) {
            typename Alloc::pointer p = get_allocator().allocate( number_depths(), number_rows()*number_columns() );
            typedef typename ecuda::add_pointer<value_type>::type raw_pointer_type;
            shared_ptr<value_type> sp( naked_cast<raw_pointer_type>(p) );
            padded_ptr< value_type, shared_ptr<value_type> > pp( sp, p.get_pitch() );
            base_type base( pp, number_rows()*number_columns(), number_depths() );
            base_type::swap( base );
        }
    }

public:
    __HOST__ cube(
        const size_type numberRows=0, const size_type numberColumns=0, const size_type numberDepths=0,
        const value_type& value = value_type(),
        const Alloc& allocator = Alloc()
    ) : base_type( pointer(), numberRows*numberColumns, numberDepths ), numberRows(numberRows), allocator(allocator)
    {
        init();
        if( size() ) ecuda::fill( begin(), end(), value );
    }

    __HOST__ cube( const cube& src ) :
        base_type( pointer(), src.number_rows()*src.number_columns(), src.number_depths() ),
        numberRows( src.numberRows ),
        allocator( src.get_allocator() )
        //allocator( std::allocator_traits<Alloc>::select_on_container_copy_construction(src.get_allocator()) )
    {
        init();
        if( size() ) ecuda::copy( src.begin(), src.end(), begin() );
    }

    __HOST__ cube( const cube& src, const allocator_type& alloc ) :
        base_type( pointer(), src.number_rows()*src.number_columns(), src.number_depths() ),
        numberRows( src.numberRows ),
        allocator(alloc)
    {
        init();
        if( size() ) ecuda::copy( src.begin(), src.end(), begin() );
    }

    __HOST__ cube& operator=( const cube& src )
    {
        if( number_rows() != src.number_rows() || number_columns() != src.number_columns() || number_depths() != src.number_depths() ) {
            resize( src.number_rows(), src.number_columns(), src.number_depths() );
        }
        if( size() ) ecuda::copy( src.begin(), src.end(), begin() );
        return *this;
    }

    #ifdef ECUDA_CPP11_AVAILABLE
    __HOST__ cube( cube&& src ) : base_type(std::move(src)), numberRows(std::move(src.numberRows)), allocator(std::move(src.allocator)) {}

    __HOST__ cube& operator=( cube&& src ) {
        base_type::operator=(std::move(src));
        allocator = std::move(src.allocator);
        numberRows = std::move(src.numberRows);
        return *this;
    }
    #endif

    __HOST__ inline allocator_type get_allocator() const { return allocator; }

    __HOST__ __DEVICE__ inline size_type number_rows() const __NOEXCEPT__ { return numberRows; }

    __HOST__ __DEVICE__ inline size_type number_columns() const __NOEXCEPT__ { return base_type::number_rows()/numberRows; } //TODO: this costs a register in kernel code

    __HOST__ __DEVICE__ inline size_type number_depths() const __NOEXCEPT__ { return base_type::number_columns(); }

    __HOST__ __DEVICE__ inline size_type size() const __NOEXCEPT__ { return base_type::size(); }

    __HOST__ __DEVICE__ inline bool empty() const __NOEXCEPT__ { return !size(); }

    __HOST__ __DEVICE__ inline pointer data() __NOEXCEPT__ { return base_type::get_pointer(); }

    __HOST__ __DEVICE__ inline const_pointer data() const __NOEXCEPT__ { return base_type::get_pointer(); }

    __HOST__ __DEVICE__ inline iterator begin() __NOEXCEPT__ { return base_type::begin(); }

    __HOST__ __DEVICE__ inline iterator end() __NOEXCEPT__ { return base_type::end(); }

    __HOST__ __DEVICE__ inline const_iterator begin() const __NOEXCEPT__ { return base_type::begin(); }

    __HOST__ __DEVICE__ inline const_iterator end() const __NOEXCEPT__ { return base_type::end(); }

    __HOST__ __DEVICE__ inline reverse_iterator rbegin() __NOEXCEPT__ { return base_type::rbegin(); }

    __HOST__ __DEVICE__ inline reverse_iterator rend() __NOEXCEPT__ { return base_type::rend(); }

    __HOST__ __DEVICE__ inline const_reverse_iterator rbegin() const __NOEXCEPT__ { return base_type::rbegin(); }

    __HOST__ __DEVICE__ inline const_reverse_iterator rend() const __NOEXCEPT__ { return base_type::rend(); }

    #ifdef ECUDA_CPP11_AVAILABLE
    __HOST__ __DEVICE__ inline const_iterator cbegin() const    __NOEXCEPT__ { return base_type::cbegin();  }
    __HOST__ __DEVICE__ inline const_iterator cend() const      __NOEXCEPT__ { return base_type::cend();    }
    __HOST__ __DEVICE__ inline const_reverse_iterator crbegin() __NOEXCEPT__ { return base_type::crbegin(); }
    __HOST__ __DEVICE__ inline const_reverse_iterator crend()   __NOEXCEPT__ { return base_type::crend();   }
    #endif

    __HOST__ __DEVICE__ row_type get_row( const size_type columnIndex, const size_type depthIndex )
    {
        typedef typename make_unmanaged<typename base_type::pointer>::type unmanaged_pointer_type;
        unmanaged_pointer_type ptr = unmanaged_cast( base_type::get_pointer() ); // padded_ptr
        ptr.skip_bytes( columnIndex*ptr.get_pitch() ); // move pointer to correct column
        ptr += depthIndex; // move pointer to correct depth
        typename row_type::pointer ptr2( ptr.get(), number_columns()*ptr.get_pitch() );
        return row_type( ptr2, number_rows() );
    }

    __HOST__ __DEVICE__ column_type get_column( const size_type rowIndex, const size_type depthIndex )
    {
        typedef typename make_unmanaged<typename base_type::pointer>::type unmanaged_pointer_type;
        unmanaged_pointer_type ptr = unmanaged_cast( base_type::get_pointer() ); // padded_ptr
        ptr.skip_bytes( rowIndex*number_columns()*ptr.get_pitch() ); // move pointer to correct column
        ptr += depthIndex; // move pointer to correct depth
        typename column_type::pointer ptr2( ptr.get(), ptr.get_pitch() );
        return column_type( ptr2, number_columns() );
    }

    __HOST__ __DEVICE__ depth_type get_depth( const size_type rowIndex, const size_type columnIndex )
    {
        typedef typename make_unmanaged<typename base_type::pointer>::type unmanaged_pointer_type;
        unmanaged_pointer_type ptr = unmanaged_cast( base_type::get_pointer() );
        ptr += rowIndex*number_columns()*number_depths()+columnIndex*number_depths(); // move pointer to depth start
        return depth_type( naked_cast<typename ecuda::add_pointer<value_type>::type>(ptr), number_depths() );
    }

    __HOST__ __DEVICE__ const_row_type get_row( const size_type columnIndex, const size_type depthIndex ) const
    {
        typedef typename make_unmanaged_const<typename base_type::pointer>::type unmanaged_pointer_type;
        unmanaged_pointer_type ptr = unmanaged_cast( base_type::get_pointer() ); // padded_ptr
        ptr.skip_bytes( columnIndex*ptr.get_pitch() ); // move pointer to correct column
        ptr += depthIndex; // move pointer to correct depth
        typename const_row_type::pointer ptr2( ptr.get(), number_columns()*ptr.get_pitch() );
        return const_row_type( ptr2, number_rows() );
    }

    __HOST__ __DEVICE__ const_column_type get_column( const size_type rowIndex, const size_type depthIndex ) const
    {
        typedef typename make_unmanaged_const<typename base_type::pointer>::type unmanaged_pointer_type;
        unmanaged_pointer_type ptr = unmanaged_cast( base_type::get_pointer() ); // padded_ptr
        ptr.skip_bytes( rowIndex*number_columns()*ptr.get_pitch() ); // move pointer to correct column
        ptr += depthIndex; // move pointer to correct depth
        typename const_column_type::pointer ptr2( ptr.get(), ptr.get_pitch() );
        return const_column_type( ptr2, number_columns() );
    }

    __HOST__ __DEVICE__ const_depth_type get_depth( const size_type rowIndex, const size_type columnIndex ) const
    {
        typedef typename make_unmanaged_const<typename base_type::pointer>::type unmanaged_pointer_type;
        unmanaged_pointer_type ptr = unmanaged_cast( base_type::get_pointer() );
        ptr += rowIndex*number_columns()*number_depths()+columnIndex*number_depths(); // move pointer to depth start
        return const_depth_type( naked_cast<typename ecuda::add_pointer<const value_type>::type>(ptr), number_depths() );
    }

    __HOST__ __DEVICE__ slice_yz_type get_yz( const size_type rowIndex )
    {
        typedef typename make_unmanaged<typename base_type::pointer>::type unmanaged_pointer_type;
        unmanaged_pointer_type ptr = unmanaged_cast( base_type::get_pointer() );
        ptr += rowIndex*number_columns()*number_depths();
        return slice_yz_type( ptr, number_columns(), number_depths() );
    }

    __HOST__ __DEVICE__ slice_xy_type get_xy( const size_type depthIndex )
    {
        typedef typename make_unmanaged<typename base_type::pointer>::type unmanaged_pointer_type;
        unmanaged_pointer_type ptr = unmanaged_cast( base_type::get_pointer() ); // padded_ptr
        ptr += depthIndex; // move to correct depth
        typename slice_xy_type::pointer ptr2( ptr.get(), ptr.get_pitch() ); // make pointer stride over depths
        return slice_xy_type( ptr2, number_rows(), number_columns() );
    }

    __HOST__ __DEVICE__ slice_xz_type get_xz( const size_type columnIndex )
    {
        typedef typename make_unmanaged<typename base_type::pointer>::type unmanaged_pointer_type;
        unmanaged_pointer_type ptr = unmanaged_cast( base_type::get_pointer() ); // padded_ptr
        typename slice_xz_type::pointer ptr2( ptr.get(), ptr.get_pitch()*number_columns() ); // extend pitch
        ptr2 += columnIndex*number_depths(); // move to correct column
        return slice_xz_type( ptr2, number_rows(), number_depths() );
    }

    __HOST__ __DEVICE__ const_slice_yz_type get_yz( const size_type rowIndex ) const
    {
        typedef typename make_unmanaged_const<typename base_type::pointer>::type unmanaged_pointer_type;
        unmanaged_pointer_type ptr = unmanaged_cast( base_type::get_pointer() );
        ptr += rowIndex*number_columns()*number_depths();
        return const_slice_yz_type( ptr, number_columns(), number_depths() );
    }

    __HOST__ __DEVICE__ const_slice_xy_type get_xy( const size_type depthIndex ) const
    {
        typedef typename make_unmanaged_const<typename base_type::pointer>::type unmanaged_pointer_type;
        unmanaged_pointer_type ptr = unmanaged_cast( base_type::get_pointer() ); // padded_ptr
        ptr += depthIndex; // move to correct depth
        typename const_slice_xy_type::pointer ptr2( ptr.get(), ptr.get_pitch() ); // make pointer stride over depths
        return const_slice_xy_type( ptr2, number_rows(), number_columns() );
    }

    __HOST__ __DEVICE__ const_slice_xz_type get_xz( const size_type columnIndex ) const
    {
        typedef typename make_unmanaged_const<typename base_type::pointer>::type unmanaged_pointer_type;
        unmanaged_pointer_type ptr = unmanaged_cast( base_type::get_pointer() ); // padded_ptr
        typename const_slice_xz_type::pointer ptr2( ptr.get(), ptr.get_pitch()*number_columns() ); // extend pitch
        ptr2 += columnIndex*number_depths(); // move to correct column
        return const_slice_xz_type( ptr2, number_rows(), number_depths() );
    }

    __DEVICE__ inline reference at( size_type rowIndex, size_type columnIndex, size_type depthIndex )
    {
        if( rowIndex >= number_rows() || columnIndex >= number_columns() || depthIndex >= number_depths() ) {
            #ifndef __CUDACC__
            throw std::out_of_range( EXCEPTION_MSG("ecuda::cube::at() row, column and/or depth index parameter is out of range") );
            #else
            // this strategy is taken from:
            // http://stackoverflow.com/questions/12521721/crashing-a-kernel-gracefully
            ecuda::threadfence();
            //__threadfence();
            asm("trap;");
            #endif
        }
        return base_type::at( rowIndex*number_columns()+columnIndex, depthIndex );
    }

    __DEVICE__ inline const_reference at( size_type rowIndex, size_type columnIndex, size_type depthIndex ) const
    {
        if( rowIndex >= number_rows() || columnIndex >= number_columns() || depthIndex >= number_depths() ) {
            #ifndef __CUDACC__
            throw std::out_of_range( EXCEPTION_MSG("ecuda::cube::at() row, column and/or depth index parameter is out of range") );
            #else
            // this strategy is taken from:
            // http://stackoverflow.com/questions/12521721/crashing-a-kernel-gracefully
            ecuda::threadfence();
            //__threadfence();
            asm("trap;");
            #endif
        }
        return base_type::at( rowIndex*number_columns()+columnIndex, depthIndex );
    }

    __DEVICE__ inline reference operator()( const size_type rowIndex, const size_type columnIndex, const size_type depthIndex ) { return base_type::at( rowIndex*number_columns()+columnIndex, depthIndex ); }

    __DEVICE__ inline const_reference operator()( const size_type rowIndex, const size_type columnIndex, const size_type depthIndex ) const { return base_type::at( rowIndex*number_columns()+columnIndex, depthIndex ); }

    __HOST__ __DEVICE__ inline slice_yz_type operator[]( const size_type rowIndex ) { return get_yz( rowIndex ); }

    __HOST__ __DEVICE__ inline const_slice_yz_type operator[]( const size_type rowIndex ) const { return get_yz( rowIndex ); }

    __HOST__ void resize( const size_type newNumberRows, const size_type newNumberColumns, const size_type newNumberDepths, const value_type& value = value_type() )
    {
        if( number_rows() == newNumberRows && number_columns() == newNumberColumns && number_depths() == newNumberDepths ) return; // no resize needed
        cube newCube( newNumberRows, newNumberColumns, newNumberDepths, value, get_allocator() );
        for( size_type i = 0; i < std::min(newNumberRows,number_rows()); ++i ) {
            for( size_type j = 0; j < std::min(newNumberColumns,number_columns()); ++j ) {
                depth_type inputDepth = get_depth(i,j);
                depth_type outputDepth = newCube.get_depth(i,j);
                ecuda::copy( inputDepth.begin(), inputDepth.begin()+std::min(newNumberDepths,number_depths()), outputDepth.begin() );
            }
        }
        base_type::swap( newCube );
        numberRows = newNumberRows;
    }

    __HOST__ __DEVICE__ inline void fill( const value_type& value ) { if( !empty() ) ecuda::fill( begin(), end(), value ); }

};

namespace impl {

template< typename T, class Alloc=device_pitch_allocator<T> >
class cube_kernel_argument : public cube<T,Alloc,typename ecuda::add_pointer<T>::type>
{

private:
    typedef cube<T,Alloc,typename ecuda::add_pointer<T>::type> base_type;

public:
    template<class P>
    __HOST__ cube_kernel_argument( const cube<T,Alloc,P>& src ) : base_type( src, ecuda::true_type() ) {}

    __HOST__ __DEVICE__ cube_kernel_argument( const cube_kernel_argument& src ) : base_type( src, ecuda::true_type() ) {}

    template<class P>
    __HOST__ cube_kernel_argument& operator=( const cube<T,Alloc,P>& src )
    {
        base_type::shallow_assign( src );
        return *this;
    }

    #ifdef ECUDA_CPP11_AVAILABLE
    cube_kernel_argument( cube_kernel_argument&& src ) : base_type(std::move(src)) {}

    cube_kernel_argument& operator=( cube_kernel_argument&& src )
    {
        base_type::operator=(std::move(src));
        return *this;
    }
    #endif

};

} // namespace impl

} // namespace ecuda

#endif