Proto_LevelBoxData.H

Go to the documentation of this file.

#pragma once
#ifndef _PROTO_LEVELBOXDATA_H_
#define _PROTO_LEVELBOXDATA_H_

#include "Proto_MemType.H"
#include "Proto_Memory.H"
#include "Proto_Point.H"
#include "Proto_Centering.H"
#include "Proto_Box.H"
#include "Proto_Centering.H"
#include "Proto_Copier.H"
#include "Proto_DisjointBoxLayout.H"
#include "Proto_DataIterator.H"
#include "Proto_NeighborIterator.H"
#include "Proto_MayDay.H"
#include "Proto_SPMD.H"
#include <cstdlib> //for size_t
#include <iostream>
#include <cstring> // for Writing data to .vtk files.

namespace Proto 
{
    typedef MotionItem<BoxPartition, BoxPartition> LevelMotionItem;
    // Forward Declarations
    template<typename T, unsigned int C, MemType MEM, Centering CTR>
    class LevelBoxData;


// =======================================================================
// LEVEL BOX DATA
    
    template<typename T, unsigned int C, MemType MEM, Centering CTR>
    class LevelExchangeCopier;

    /// Level Box Data
    /**
        A distributed data holder implementing the abstraction of a union of logically rectangular patches.
    */
    template
    <typename T, unsigned int C=1, MemType MEM=MEMTYPE_DEFAULT, Centering CTR = PR_CELL>
    class LevelBoxData
    {
        public:

        inline LevelBoxData();

        inline LevelBoxData(const DisjointBoxLayout& a_layout, const Point& a_ghostRegionSizes);

        template<unsigned int DST_DATA_COMPONENTS>
        inline LevelBoxData(const LevelBoxData<T, DST_DATA_COMPONENTS, MEM, CTR>& a_srcData, unsigned int a_firstComponent);

        inline void define(const DisjointBoxLayout& a_layout, const Point& a_ghostRegionSizes);

        inline unsigned int numBoxes() const;

        inline void setToZero(int a_comp = -1);
        
        inline void setVal(T a_value, int a_comp = -1);

        /** Set data to random noise between two values */
        inline void setRandom(T a_low, T a_high); 

        /** Periodic boundaries are not considered domain boundaries */
        inline void setDomainBoundary(T a_value, int a_comp = -1);
      
        inline void iota(Array<T,DIM> a_dx, Array<T,DIM> a_offset = Array<T,DIM>::Zeros());

        /// Function Initialization
        /**
            Initializes *this using a function with a signature similar to that used by
            Proto::forallInPlace_p. 

            Func signature:
            void a_func(Point& a_point, Var<T, C, MEM>& a_thisData, ... )

            Like Proto::forall and its variants, this function can accept additional data
            holders (e.g. other LevelBoxData instances) so long as all inputs are compatible
            in the sense of LevelBoxData::compatible(...).
        */
        template<typename Func, typename... Srcs>
        inline void initialize(Func& a_func, Srcs&... a_srcs);
        
        /// Exchange
        /**
            Copies data from the valid regions of *this into ghost regions. When MPI is
            enabled, this function also takes care of any necessary communication between
            patches on different processes.
        */
        inline void exchange();

        template< template<typename, unsigned int, MemType, Centering> class E_COPIER,
            typename... Args>
        inline void defineExchange(Args... a_args);
        /// Copy To
        /**
            Copies data from the valid regions of this into the valid regions of another
            LevelBoxData. This function copies both valid data and ghost regions if applicable.
            If <code>a_dest</code> has a larger ghost region then <code>*this</code>, the user
            is expected to call <code>exchange</code> to correctly fill the extra ghost cells.

            \param a_dest   The destination to which data in *this will be copied.
        */
        template<MemType DST_MEM>
        inline void copyTo(LevelBoxData<T, C, DST_MEM, CTR>& a_dest) const;
        
        template<MemType DST_MEM>
        inline void copyToSimple(LevelBoxData<T, C, DST_MEM, CTR>& a_dest) const;
       
        /// Copy With Coarsen
        /**
            Coarsen's the data in *this by applying a geometric average before copying 
            into the destination.
        */
        inline void coarsenTo(LevelBoxData& a_dest, Point a_refRatio) const;

        /// Linear Offset
        /**
          Returns the serial index of the first data element of this stored on
          this a_proc.
         */
        inline unsigned int offset(int a_proc) const;
        inline unsigned int offset() const { return offset(Proto::procID()); }

        /// Patch Size
        /**
          Computes the size (in data cells) of a single patch in *this.
          The output accounts for both ghost cells and centering.
         */
        inline unsigned int patchSize() const;
        
        /// Patch Box
        /**
            Returns the Box corresponding to the valid cells of the data at a given index.
            This box does not include ghost cells, but does account for centering.
            For cell centered data, this is the same as layout()[a_index], but will
            in general be different for other types of centering.
         */
        inline Box patchBox(const DataIndex<BoxPartition>& a_index) const;

        /// Linear In (Patch, Box)
        /**
          Reads data from a serial buffer, populating the patch associated with 
          a_index within the subset specified by a_box. 

          \param a_buffer     The read buffer
          \param a_box        A box defining the subset to read into
          \param a_index      The index defining the patch to read into
         */ 
        inline void linearIn(void* a_buffer, 
                const Box& a_box,
                const DataIndex<BoxPartition>& a_index);

        /// Linear Out (Patch, Box)
        /**
          Writes data to a serial buffer, from the patch associated with 
          a_index within the subset specified by a_box. 

          \param a_buffer     The write buffer
          \param a_box        A box defining the subset to write from
          \param a_index      The index defining the patch to write from
         */ 
        inline void linearOut(void* a_buffer, 
                const Box& a_box,
                const DataIndex<BoxPartition>& a_index) const;

        /// Serial Size (Patch, Box)
        /**
          Computes the size in bytes of the serial data buffer needed to hold
          the patch associated with a_index within the subset specified by a_box.

          \param a_box        A box defining the subset of a patch
          \param a_index      The index defining the patch
         */ 
        inline size_t linearSize(const Box& a_box,
                const DataIndex<BoxPartition>& a_index) const;
        /// Linear In (Patch)
        /**
          Overload of LinearIn that always reads in a full patch
         */
        inline void linearIn(void* a_buffer, 
                const DataIndex<BoxPartition>& a_index);
        /// Linear Out (Patch)
        /**
          Overload of LinearOut that always writes out a full patch
         */
        inline void linearOut(void* a_buffer, 
                const DataIndex<BoxPartition>& a_index) const;
        /// Size (Patch)
        /**
          Overload of Size that always computes the size of a full patch
         */
        inline size_t linearSize( const DataIndex<BoxPartition>& a_index) const;

        /// Linear In (All Local Data)
        /**
          Overload of LinearIn that reads all local data
         */
        inline void linearIn(void* a_buffer);
        /// Linear Out (All Local Data)
        /**
          Overload of LinearOut that writes out all local data
         */
        inline void linearOut(void* a_buffer) const;
        /// Size (All Local Data)
        /**
          Overload of Size that computes the full local size of the data in *this
         */
        inline size_t linearSize() const;

        /// Get Box Layout
        inline DisjointBoxLayout layout() const { return m_layout; }

        /// Get Ghost Size
        inline Point ghost() const { return m_ghost; }
       
        /// Increment
        /**
            Add a scaled multiple of another LevelBoxData.
        */
        inline void increment(LevelBoxData<T, C, MEM, CTR>& a_data, T a_scale = 1.0);
       
        /// Multiply By Scalar
        inline void operator*=(T a_scale);

        /// Add a Constant
        inline void operator+=(T a_scale);
        
        /// Reduction
        /**
            Computes a reduction for a component over all valid cells (ghost region is ignored)
        */
        template<Proto::Operation OP>
        inline double reduce(unsigned int a_comp = 0) const;
       
        /// Maximum Absolute Value
        /**
            Computes the maximum absolute value for a component over all valid cells (ghost region is ignored)
        */
        inline double absMax(unsigned int a_comp = 0) const;
        
        /// Sum
        /**
            Computes the sum of a component over all valid cells (ghost region is ignored)
        */
        inline double sum(unsigned int a_comp = 0) const;
        
        /// Max
        /**
            Computes the max value of a component over all valid cells (ghost region is ignored)
        */
        inline double max(unsigned int a_comp = 0) const;
        
        /// Min
        /**
            Computes the min value of a component over all valid cells (ghost region is ignored)
        */
        inline double min(unsigned int a_comp = 0) const;
        
        /// Integral
        /**
            Compute the integral over this using an isotropic grid spacing.

            \param a_dx     Isotropic grid spacing
            \param a_comp   A component to integrate over
        */
        inline double integrate(double a_dx, unsigned int a_comp = 0) const;
        
        /// Integral (Anisotropic)
        /**
            Compute the integral over this using an isotropic grid spacing.

            \param a_dx     Isotropic grid spacing
            \param a_comp   A component to integrate over
        */
        inline double integrate(Array<double, DIM> a_dx, unsigned int a_comp = 0) const;
        
        /// Integral of Absolute Value
        /**
            Compute the integral of the absolute value over this using an isotropic grid spacing.

            \param a_dx     Isotropic grid spacing
            \param a_comp   A component to integrate over
        */
        inline double integrateAbs(double a_dx, unsigned int a_comp = 0) const;
        
        /// Integral of Absolute Value (Anisotropic)
        /**
            Compute the integral of the absolute value over this using an isotropic grid spacing.

            \param a_dx     Isotropic grid spacing
            \param a_comp   A component to integrate over
        */
        inline double integrateAbs(Array<double, DIM> a_dx, unsigned int a_comp = 0) const;
      
        /// Query Layout Compatibility
        template<typename _T, unsigned int _C, MemType _MEM, Centering _CTR>
        inline bool compatible(const LevelBoxData<_T, _C, _MEM, _CTR>& a_data) const;

        /// Query Layout Compatibility
        inline bool compatible(const DisjointBoxLayout& a_layout) const;
        /// Query Layout Compatibility
        inline bool compatible(const DataIterator<BoxPartition>& a_iter) const;
        /// Query Layout Compatibility
        inline bool compatible(const DataIndex<BoxPartition>& a_index) const;

        //=======================================
        // DATA ITERATOR INTERFACE

        /// Get Patch
        inline BoxData<T,C,MEM>& operator[](const DataIndex<BoxPartition>& a_index);
        /// Get Const Patch
        inline const BoxData<T,C,MEM>& operator[](const DataIndex<BoxPartition>& a_index) const;
        /// Get Patch (Face / Edge Centering)
        inline BoxData<T,C,MEM>& operator()(const DataIndex<BoxPartition>& a_index, unsigned int a_dir = 0);
        /// Get Const Patch (Face / Edge Centering)
        inline const BoxData<T,C,MEM>& operator()(const DataIndex<BoxPartition>& a_index, unsigned int a_dir = 0) const;
        
        /// Get Iterator
        inline DataIterator<BoxPartition> begin() const;
        inline DataIterator<BoxPartition> end() const;

      /// deprecated  but useful
      inline DataIterator<BoxPartition> dataIterator() const
      {
        return begin();
      }
        /// For debugging purposes.
        int s_verbosity = 0;
        
        private: 

        std::vector<std::vector<shared_ptr<BoxData<T, C, MEM> > >> m_data;
        Point                                   m_ghost;
        DisjointBoxLayout                       m_layout;
        bool                                    m_isDefined;
        std::shared_ptr<LevelExchangeCopier<T, C, MEM, CTR>>     m_exchangeCopier;
        
        // disallow copy constructors and assignment operators
        // to avoid very hard-to-find performance problems
        LevelBoxData<T, C, MEM, CTR>& operator= (const LevelBoxData<T, C, MEM, CTR>& a_rhs);
        LevelBoxData<T, C, MEM, CTR>(const LevelBoxData<T, C, MEM, CTR>& a_rhs);
    }; // end class LevelBoxData
    
// =======================================================================
// LEVEL BOX DATA COPIER OP
    
    /// Level Copier Operator
    /**
        Copier operator which provides the necessary interface to copy to and from
        LevelBoxData using the Copier interface. See Proto_Copier.H.
    */
    template
    <typename T, unsigned int C, MemType SRC_MEM, MemType DST_MEM, Centering CTR>
    class LevelCopierOp
    {
        public:
        LevelBoxData<T,C,SRC_MEM,CTR>*  m_src;
        LevelBoxData<T,C,DST_MEM,CTR>*  m_dst;

        inline LevelCopierOp(){};
        
        inline LevelCopierOp(
            LevelBoxData<T, C, SRC_MEM, CTR>& a_src,
            LevelBoxData<T, C, DST_MEM, CTR>& a_dst);
        
        /// Get Serialized Size
        inline int linearSize(const Box& a_bx,
                const DataIndex<BoxPartition>& a_index) const;
        
        /// Serial Write To Buffer
        inline void linearOut(void* a_buf, const LevelMotionItem& a_info);
        //inline void linearOut(void* a_buf, 
        //        const Box& a_bx,
        //        const DataIndex<BoxPartition>& a_index) const;
        
        /// Serial Read From Buffer
        inline void linearIn(void* a_buf, const LevelMotionItem& a_info);
        //inline void linearIn(void* a_buf, 
        //        const Box& a_bx,
        //        const DataIndex<BoxPartition>& a_index) const;
        
        /// Local Copy
        inline void localCopy(const LevelMotionItem& a_info);
        //inline void localCopy(
        //    const Box& a_domain, const DataIndex<BoxPartition>& a_domainIndex,
        //    const Box& a_range,  const DataIndex<BoxPartition>& a_rangeIndex) const;
    };
// =======================================================================
// LEVEL BOX DATA COPIERS
    
    /// Level Copier
    /**
        A Copier used to implement copying of data between LevelBoxData objects.
    */
    template<typename T, unsigned int C, MemType SRC_MEM, MemType DST_MEM, Centering CTR>
    class LevelCopier 
    : public Copier<LevelCopierOp<T, C, SRC_MEM, DST_MEM, CTR>, BoxPartition, BoxPartition, SRC_MEM, DST_MEM>
    {
        public:
        
        inline LevelCopier(){};
        
        /// Build Copier Motion Plan
        inline void buildMotionPlans(LevelCopierOp<T, C, SRC_MEM, DST_MEM, CTR>& a_op);
    };
    
    /// Exchange Copier
    /**
        A Copier used to execute the data copies necessary for filling ghost regions
        within a LevelBoxData.
    */
    template<typename T, unsigned int C, MemType MEM, Centering CTR>
    class LevelExchangeCopier
    : public Copier<LevelCopierOp<T, C, MEM, MEM, CTR>, BoxPartition, BoxPartition, MEM, MEM>
    {
        public:

        inline LevelExchangeCopier(){};
        
        /// Build Copier Motion Plan
        inline virtual void buildMotionPlans(LevelCopierOp<T, C, MEM, MEM, CTR>& a_op);
    };


    /// Interpolate Boundaries
    /**
        General utility function for interpolating data into coarse-fine boundaries.
        This version creates a temporary LevelBoxData representing the coarsened fine region.
    */
    template<typename T, unsigned int C, MemType MEM, Centering CTR>
    void interpBoundaries(
        LevelBoxData<T, C, MEM, CTR>& a_crse,
        LevelBoxData<T, C, MEM, CTR>& a_fine,
        InterpStencil<T>&       a_interp);
    
    /// Interpolate Boundaries
    /**
        General utility function for interpolating data into coarse-fine boundaries.
        This version does not create any temporary LevelBoxData but requires an additional
        input dataholder representing the coarsened fine region.
    */
    template<typename T, unsigned int C, MemType MEM, Centering CTR>
    void interpBoundaries(
        LevelBoxData<T, C, MEM, CTR>& a_crse,
        LevelBoxData<T, C, MEM, CTR>& a_fine,
        LevelBoxData<T, C, MEM, CTR>& a_crseFine,
        InterpStencil<T>&       a_interp);

    /// Average Down
    /**
        General utility function for averaging fine data onto coarse data
    */
    template<typename T, unsigned int C, MemType MEM, Centering CTR>
    void averageDown(
        LevelBoxData<T, C, MEM, CTR>& a_crse,
        LevelBoxData<T, C, MEM, CTR>& a_fine,
        Point                         a_refRatio);
    
    /// Average Down (Scalar Refinement Ratio)
    /**
        General utility function for averaging fine data onto coarse data
    */
    template<typename T, unsigned int C, MemType MEM, Centering CTR>
    void averageDown(
        LevelBoxData<T, C, MEM, CTR>& a_crse,
        LevelBoxData<T, C, MEM, CTR>& a_fine,
        int                           a_refRatio);
    
    /// Average Down
    /**
        General utility function for averaging fine data onto coarse data
    */
    template<typename T, unsigned int C, MemType MEM, Centering CTR>
    void averageDown(
        LevelBoxData<T, C, MEM, CTR>& a_crse,
        LevelBoxData<T, C, MEM, CTR>& a_fine,
        LevelBoxData<T, C, MEM, CTR>& a_crseFine,
        Point                         a_refRatio);
    
    /// Average Down (Scalar Refinement Ratio)
    /**
        General utility function for averaging fine data onto coarse data
    */
    template<typename T, unsigned int C, MemType MEM, Centering CTR>
    void averageDown(
        LevelBoxData<T, C, MEM, CTR>& a_crse,
        LevelBoxData<T, C, MEM, CTR>& a_fine,
        LevelBoxData<T, C, MEM, CTR>& a_crseFine,
        int                           a_refRatio);
    

#include "implem/Proto_LevelBoxDataImplem.H"
} //end namespace Proto
#endif