/**
 * @file   grid.hpp
 * @author Julian Iseringhausen <isering@ins.uni-bonn.de>
 * @date   Mon Apr 18 12:53:45 2011
 *
 * @brief  VMG::Grid
 *
 */

#ifndef GRID_HPP_
#define GRID_HPP_

#include "base/object.hpp"
#include "grid/grid_indexing.hpp"

namespace VMG
{

class Comm;

class Grid : public Object
{
public:
  Grid()
  {
    grid = NULL;
  }

  Grid(const GlobalIndices& global_, const LocalIndices& local_, const SpatialExtent& extent_) :
    global(global_),
    local(local_),
    extent(extent_)
  {InitGrid();}

  virtual ~Grid();

  Grid(const Comm* comm, const Index& size, const vmg_float& width);

  const GlobalIndices& Global() const {return global;}
  const LocalIndices& Local() const {return local;}
  const SpatialExtent& Extent() const {return extent;}

  GlobalIndices& Global() {return global;}
  LocalIndices& Local() {return local;}
  SpatialExtent& Extent() {return extent;}

  const vmg_float& MeshWidth() const {return Extent().MeshWidth().X();} ///< Mesh width of current level

  static vmg_float& Correction() {return Grid::correction;}

  void Clear();         ///< Overwrites all grid points on current level with zeros
  void ClearInner();
  void ClearHalo();     ///< Overwrites all halo points on current level with zeros
  void ClearBoundary(); ///< Overwrites all boundary points on current level with zeros

  vmg_float& operator()(int x, int y, int z);  ///< Returns a reference to the requested gridpoint.
  vmg_float& operator()(const Index& index);

  const vmg_float& GetVal(int x, int y, int z) const; ///< Returns the value of a requested gridpoint.
  const vmg_float& GetVal(const Index& index) const;

  vmg_float GetCorrectedVal(int x, int y, int z) const {return this->GetVal(x, y, z) - correction;} ///< Return grid
  vmg_float GetCorrectedVal(const Index& index) const {return this->GetVal(index) - correction;}

  void ForceDiscreteCompatibilityCondition();
  void SetAverageToZero();

  void SetGrid(const Grid& rhs);     ///< Overwrite current grid with values from another grid
  void SetBoundary(const Grid& rhs); ///< Overwrite boundary with values from rhs

  void AddGrid(const Grid& rhs);             ///< Add values of another grid
  void SubtractGrid(const Grid& rhs);        ///< Subtract values of another grid
  void MultiplyScalar(const vmg_float& scalar);   ///< Multiply grid values with scalar
  void ApplyStencil(const Stencil& stencil); ///< Apply stencil to grid

  int GlobalLinearIndex(int x, int y, int z) const;  ///< Compute a unique 1-dimensional global index
  int GlobalLinearIndex(const Index& index) const;

  Index GlobalToLocalIndex(const Index& index) const;

  bool IsCompatible(const Grid& rhs) const; ///< Check if two grids share compatible settings
  bool IsConsistent() const;                  ///< Check grid for nan and inf

  virtual vmg_float DebugKnownSolution(Vector x) const {return 0.0;}
  vmg_float DebugKnownSolution(Index i) const {return DebugKnownSolution(extent.Begin() + i * extent.MeshWidth());}

private:
  void InitGrid();

protected:
  GlobalIndices global;
  LocalIndices local;
  SpatialExtent extent;

  vmg_float *grid;

  static vmg_float correction;
};

inline vmg_float& Grid::operator()(int x, int y, int z)
{
  return grid[z + local.SizeTotal().Z() * (y + local.SizeTotal().Y() * x)];
}

inline vmg_float& Grid::operator()(const Index& index)
{
  return this->operator()(index.X(), index.Y(), index.Z());
}

inline const vmg_float& Grid::GetVal(int x, int y, int z) const
{
  return grid[z + local.SizeTotal().Z() * (y + local.SizeTotal().Y() * x)];
}

inline const vmg_float& Grid::GetVal(const Index& index) const
{
  return this->GetVal(index.X(), index.Y(), index.Z());
}

inline int Grid::GlobalLinearIndex(int x, int y, int z) const
{
  return z + global.Size().Z() * (y + global.Size().Y() * x);
}

inline int Grid::GlobalLinearIndex(const Index& index) const
{
  return GlobalLinearIndex(index.X(), index.Y(), index.Z());
}

inline Index Grid::GlobalToLocalIndex(const Index& index) const
{
  return index - Global().Begin() + Local().Begin();
}

}

#endif /* GRID_HPP_ */