| [de061d] | 1 | /*
 | 
|---|
 | 2 |  *    vmg - a versatile multigrid solver
 | 
|---|
 | 3 |  *    Copyright (C) 2012 Institute for Numerical Simulation, University of Bonn
 | 
|---|
 | 4 |  *
 | 
|---|
 | 5 |  *  vmg is free software: you can redistribute it and/or modify
 | 
|---|
 | 6 |  *  it under the terms of the GNU General Public License as published by
 | 
|---|
 | 7 |  *  the Free Software Foundation, either version 3 of the License, or
 | 
|---|
 | 8 |  *  (at your option) any later version.
 | 
|---|
 | 9 |  *
 | 
|---|
 | 10 |  *  vmg is distributed in the hope that it will be useful,
 | 
|---|
 | 11 |  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 | 
|---|
 | 12 |  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 | 
|---|
 | 13 |  *  GNU General Public License for more details.
 | 
|---|
 | 14 |  *
 | 
|---|
 | 15 |  *  You should have received a copy of the GNU General Public License
 | 
|---|
 | 16 |  *  along with this program.  If not, see <http://www.gnu.org/licenses/>.
 | 
|---|
 | 17 |  */
 | 
|---|
 | 18 | 
 | 
|---|
 | 19 | /**
 | 
|---|
 | 20 |  * @file   grid.hpp
 | 
|---|
 | 21 |  * @author Julian Iseringhausen <isering@ins.uni-bonn.de>
 | 
|---|
 | 22 |  * @date   Mon Apr 18 12:53:45 2011
 | 
|---|
 | 23 |  *
 | 
|---|
 | 24 |  * @brief  VMG::Grid
 | 
|---|
 | 25 |  *
 | 
|---|
 | 26 |  */
 | 
|---|
 | 27 | 
 | 
|---|
 | 28 | #ifndef GRID_HPP_
 | 
|---|
 | 29 | #define GRID_HPP_
 | 
|---|
 | 30 | 
 | 
|---|
 | 31 | #include "base/object.hpp"
 | 
|---|
 | 32 | #include "grid/grid_iterator.hpp"
 | 
|---|
 | 33 | #include "grid/grid_iterator_suite.hpp"
 | 
|---|
 | 34 | #include "grid/grid_properties.hpp"
 | 
|---|
 | 35 | 
 | 
|---|
 | 36 | namespace VMG
 | 
|---|
 | 37 | {
 | 
|---|
 | 38 | 
 | 
|---|
 | 39 | class Comm;
 | 
|---|
 | 40 | class Multigrid;
 | 
|---|
 | 41 | class Stencil;
 | 
|---|
 | 42 | 
 | 
|---|
 | 43 | class Grid : public Object
 | 
|---|
 | 44 | {
 | 
|---|
 | 45 | public:
 | 
|---|
 | 46 |   typedef GridIterator iterator;
 | 
|---|
 | 47 | 
 | 
|---|
 | 48 |   Grid(int level_ = 0, Multigrid* father_ = NULL) :
 | 
|---|
 | 49 |     level(level_),
 | 
|---|
 | 50 |     father(father_)
 | 
|---|
 | 51 |   {
 | 
|---|
 | 52 |     grid = NULL;
 | 
|---|
 | 53 |   }
 | 
|---|
 | 54 | 
 | 
|---|
 | 55 |   Grid(const GlobalIndices& global_,
 | 
|---|
 | 56 |        const LocalIndices& local_,
 | 
|---|
 | 57 |        const SpatialExtent& extent_,
 | 
|---|
 | 58 |        int level_ = 0,
 | 
|---|
 | 59 |        Multigrid* father_ = NULL) :
 | 
|---|
 | 60 |     level(level_),
 | 
|---|
 | 61 |     global(global_),
 | 
|---|
 | 62 |     local(local_),
 | 
|---|
 | 63 |     extent(extent_),
 | 
|---|
 | 64 |     iterators(local_),
 | 
|---|
 | 65 |     father(father_)
 | 
|---|
 | 66 |   {
 | 
|---|
 | 67 |     InitGrid();
 | 
|---|
 | 68 |   }
 | 
|---|
 | 69 | 
 | 
|---|
 | 70 |   Grid(const Grid& rhs) :
 | 
|---|
 | 71 |     level(rhs.Level()),
 | 
|---|
 | 72 |     global(rhs.Global()),
 | 
|---|
 | 73 |     local(rhs.Local()),
 | 
|---|
 | 74 |     extent(rhs.Extent()),
 | 
|---|
 | 75 |     iterators(rhs.Iterators()),
 | 
|---|
 | 76 |     father(rhs.Father())
 | 
|---|
 | 77 |   {
 | 
|---|
 | 78 |     InitGrid();
 | 
|---|
 | 79 |     SetGrid(rhs);
 | 
|---|
 | 80 |   }
 | 
|---|
 | 81 | 
 | 
|---|
 | 82 |   virtual ~Grid();
 | 
|---|
 | 83 | 
 | 
|---|
 | 84 |   Grid& operator=(const Grid& rhs);
 | 
|---|
 | 85 | 
 | 
|---|
 | 86 |   GlobalIndices& Global() {return global;}
 | 
|---|
 | 87 |   LocalIndices& Local() {return local;}
 | 
|---|
 | 88 |   SpatialExtent& Extent() {return extent;}
 | 
|---|
 | 89 | 
 | 
|---|
 | 90 |   const GlobalIndices& Global() const {return global;}
 | 
|---|
 | 91 |   const LocalIndices& Local() const {return local;}
 | 
|---|
 | 92 |   const SpatialExtent& Extent() const {return extent;}
 | 
|---|
 | 93 | 
 | 
|---|
 | 94 |   GridIteratorSuite& Iterators() {return iterators;}
 | 
|---|
 | 95 |   const GridIteratorSuite& Iterators() const {return iterators;}
 | 
|---|
 | 96 | 
 | 
|---|
 | 97 |   void Clear();         ///< Overwrites all grid points on current level with zeros
 | 
|---|
 | 98 |   void ClearInner();
 | 
|---|
 | 99 |   void ClearHalo();     ///< Overwrites all halo points on current level with zeros
 | 
|---|
 | 100 |   void ClearBoundary(); ///< Overwrites all boundary points on current level with zeros
 | 
|---|
 | 101 | 
 | 
|---|
 | 102 |   vmg_float& operator()(int x, int y, int z);  ///< Returns a reference to the requested gridpoint.
 | 
|---|
 | 103 |   vmg_float& operator()(const Index& index);
 | 
|---|
 | 104 | 
 | 
|---|
 | 105 |   const vmg_float& operator()(int x, int y, int z) const;  ///< Returns a reference to the requested gridpoint.
 | 
|---|
 | 106 |   const vmg_float& operator()(const Index& index) const;
 | 
|---|
 | 107 | 
 | 
|---|
 | 108 |   const vmg_float& GetVal(int x, int y, int z) const; ///< Returns the value of a requested gridpoint.
 | 
|---|
 | 109 |   const vmg_float& GetVal(const Index& index) const;
 | 
|---|
 | 110 | 
 | 
|---|
 | 111 |   void ForceDiscreteCompatibilityCondition();
 | 
|---|
 | 112 |   void SetAverageToZero();
 | 
|---|
 | 113 | 
 | 
|---|
 | 114 |   void SetGrid(const Grid& rhs);     ///< Overwrite current grid with values from another grid
 | 
|---|
 | 115 |   void SetBoundary(const Grid& rhs); ///< Overwrite boundary with values from rhs
 | 
|---|
 | 116 | 
 | 
|---|
 | 117 |   void AddGrid(const Grid& rhs);             ///< Add values of another grid
 | 
|---|
 | 118 |   void SubtractGrid(const Grid& rhs);        ///< Subtract values of another grid
 | 
|---|
 | 119 |   void MultiplyScalar(const vmg_float& scalar);   ///< Multiply grid values with scalar
 | 
|---|
 | 120 |   void ApplyStencil(const Stencil& stencil); ///< Apply stencil to grid
 | 
|---|
 | 121 | 
 | 
|---|
 | 122 |   int GlobalLinearIndex(const int& x, const int& y, const int& z) const;  ///< Compute a unique 1-dimensional global index
 | 
|---|
 | 123 |   int GlobalLinearIndex(const Index& index) const;
 | 
|---|
 | 124 | 
 | 
|---|
 | 125 |   bool IsCompatible(const Grid& rhs) const; ///< Check if two grids share compatible settings
 | 
|---|
 | 126 |   bool IsConsistent() const;                  ///< Check grid for nan and inf
 | 
|---|
 | 127 | 
 | 
|---|
 | 128 |   Multigrid* Father() const {return father;}
 | 
|---|
 | 129 | 
 | 
|---|
 | 130 |   const int& Level() const {return level;}
 | 
|---|
 | 131 | 
 | 
|---|
 | 132 |   bool IsActive() const {return Local().Size().Product() > 0;}
 | 
|---|
 | 133 | 
 | 
|---|
 | 134 |   Vector GetSpatialPos(const Index& index_local) const;
 | 
|---|
 | 135 |   Vector GetSpatialPosGlobal(const Index& index_global) const;
 | 
|---|
 | 136 |   Index GetGlobalIndex(const Vector& pos) const;
 | 
|---|
 | 137 | 
 | 
|---|
 | 138 | private:
 | 
|---|
 | 139 |   void InitGrid();
 | 
|---|
 | 140 | 
 | 
|---|
 | 141 | protected:
 | 
|---|
 | 142 |   int level;
 | 
|---|
 | 143 | 
 | 
|---|
 | 144 |   GlobalIndices global;
 | 
|---|
 | 145 |   LocalIndices local;
 | 
|---|
 | 146 |   SpatialExtent extent;
 | 
|---|
 | 147 | 
 | 
|---|
 | 148 |   GridIteratorSuite iterators;
 | 
|---|
 | 149 | 
 | 
|---|
 | 150 |   vmg_float *grid;
 | 
|---|
 | 151 | 
 | 
|---|
 | 152 |   Multigrid* father;
 | 
|---|
 | 153 | };
 | 
|---|
 | 154 | 
 | 
|---|
 | 155 | inline vmg_float& Grid::operator()(int x, int y, int z)
 | 
|---|
 | 156 | {
 | 
|---|
 | 157 |   return grid[z + local.SizeTotal().Z() * (y + local.SizeTotal().Y() * x)];
 | 
|---|
 | 158 | }
 | 
|---|
 | 159 | 
 | 
|---|
 | 160 | inline vmg_float& Grid::operator()(const Index& index)
 | 
|---|
 | 161 | {
 | 
|---|
 | 162 |   return this->operator()(index.X(), index.Y(), index.Z());
 | 
|---|
 | 163 | }
 | 
|---|
 | 164 | 
 | 
|---|
 | 165 | inline const vmg_float& Grid::operator()(int x, int y, int z) const
 | 
|---|
 | 166 | {
 | 
|---|
 | 167 |   return grid[z + local.SizeTotal().Z() * (y + local.SizeTotal().Y() * x)];
 | 
|---|
 | 168 | }
 | 
|---|
 | 169 | 
 | 
|---|
 | 170 | inline const vmg_float& Grid::operator()(const Index& index) const
 | 
|---|
 | 171 | {
 | 
|---|
 | 172 |   return this->operator()(index.X(), index.Y(), index.Z());
 | 
|---|
 | 173 | }
 | 
|---|
 | 174 | 
 | 
|---|
 | 175 | inline const vmg_float& Grid::GetVal(int x, int y, int z) const
 | 
|---|
 | 176 | {
 | 
|---|
 | 177 |   return grid[z + local.SizeTotal().Z() * (y + local.SizeTotal().Y() * x)];
 | 
|---|
 | 178 | }
 | 
|---|
 | 179 | 
 | 
|---|
 | 180 | inline const vmg_float& Grid::GetVal(const Index& index) const
 | 
|---|
 | 181 | {
 | 
|---|
 | 182 |   return this->GetVal(index.X(), index.Y(), index.Z());
 | 
|---|
 | 183 | }
 | 
|---|
 | 184 | 
 | 
|---|
 | 185 | inline int Grid::GlobalLinearIndex(const int& x, const int& y, const int& z) const
 | 
|---|
 | 186 | {
 | 
|---|
 | 187 |   return z - global.GlobalBegin().X()
 | 
|---|
 | 188 |     + global.GlobalSize().Z() * (y - global.GlobalBegin().Y()
 | 
|---|
 | 189 |     + global.GlobalSize().Y() * (x - global.GlobalBegin().X()));
 | 
|---|
 | 190 | }
 | 
|---|
 | 191 | 
 | 
|---|
 | 192 | inline int Grid::GlobalLinearIndex(const Index& index) const
 | 
|---|
 | 193 | {
 | 
|---|
 | 194 |   return GlobalLinearIndex(index.X(), index.Y(), index.Z());
 | 
|---|
 | 195 | }
 | 
|---|
 | 196 | 
 | 
|---|
 | 197 | }
 | 
|---|
 | 198 | 
 | 
|---|
 | 199 | #endif /* GRID_HPP_ */
 | 
|---|