mantaflow  0.10
A framework for fluid simulation
grid.h
Go to the documentation of this file.
1 
2 /******************************************************************************
3  *
4  * MantaFlow fluid solver framework
5  * Copyright 2011 Tobias Pfaff, Nils Thuerey
6  *
7  * This program is free software, distributed under the terms of the
8  * GNU General Public License (GPL)
9  * http://www.gnu.org/licenses
10  *
11  * Grid representation
12  *
13  ******************************************************************************/
14 
15 #ifndef _GRID_H
16 #define _GRID_H
17 
18 #include "manta.h"
19 #include "vectorbase.h"
20 #include "interpol.h"
21 #include "interpolHigh.h"
22 #include "kernel.h"
23 
24 namespace Manta {
25 class LevelsetGrid;
26 
29 PYTHON class GridBase : public PbClass {
30 public:
31  enum GridType { TypeNone = 0, TypeReal = 1, TypeInt = 2, TypeVec3 = 4, TypeMAC = 8, TypeLevelset = 16, TypeFlags = 32 };
32 
33  PYTHON() GridBase(FluidSolver* parent);
34 
36  inline int getSizeX() const { return mSize.x; }
38  inline int getSizeY() const { return mSize.y; }
40  inline int getSizeZ() const { return mSize.z; }
42  inline Vec3i getSize() const { return mSize; }
43 
45  inline IndexInt getStrideX() const { return 1; }
47  inline IndexInt getStrideY() const { return mSize.x; }
49  inline IndexInt getStrideZ() const { return mStrideZ; }
50 
51  inline Real getDx() { return mDx; }
52 
54  inline void checkIndex(int i, int j, int k) const;
56  inline void checkIndex(IndexInt idx) const;
58  inline bool isInBounds(const Vec3i& p, int bnd) const;
60  inline bool isInBounds(const Vec3i& p) const;
62  inline bool isInBounds(const Vec3& p, int bnd = 0) const { return isInBounds(toVec3i(p), bnd); }
64  inline bool isInBounds(IndexInt idx) const;
65 
67  inline GridType getType() const { return mType; }
69  inline bool is3D() const { return m3D; }
70 
72  inline IndexInt index(int i, int j, int k) const { DEBUG_ONLY(checkIndex(i,j,k)); return (IndexInt)i + (IndexInt)mSize.x * j + (IndexInt)mStrideZ * k; }
74  inline IndexInt index(const Vec3i& pos) const { DEBUG_ONLY(checkIndex(pos.x,pos.y,pos.z)); return (IndexInt)pos.x + (IndexInt)mSize.x * pos.y + (IndexInt)mStrideZ * pos.z; }
75 
77  inline bool is4D() const { return false; }
78  inline int getSizeT() const { return 1; }
79  inline int getStrideT() const { return 0; }
80  inline int index(int i, int j, int k, int unused) const { return index(i,j,k); }
81  inline bool isInBounds(int i,int j, int k, int t, int bnd) const { if(t!=0) return false; return isInBounds( Vec3i(i,j,k), bnd ); }
82 
83 protected:
84 
85  GridType mType;
86  Vec3i mSize;
87  Real mDx;
88  bool m3D;
89  // precomputed Z shift: to ensure 2D compatibility, always use this instead of sx*sy !
90  IndexInt mStrideZ;
91 };
92 
95 PYTHON template<class T> class Grid : public GridBase {
96 public:
98  PYTHON() Grid(FluidSolver* parent, bool show = true);
100  Grid(const Grid<T>& a);
102  virtual ~Grid();
103 
104  typedef T BASETYPE;
105  typedef GridBase BASETYPE_GRID;
106 
107  PYTHON() void save(std::string name);
108  PYTHON() void load(std::string name);
109 
111  PYTHON() void clear();
112 
115  inline T get(int i,int j, int k) const { return mData[index(i,j,k)]; }
117  inline T& get(int i,int j, int k) { return mData[index(i,j,k)]; }
119  inline T get(IndexInt idx) const { DEBUG_ONLY(checkIndex(idx)); return mData[idx]; }
121  inline T get(const Vec3i& pos) const { return mData[index(pos)]; }
123  inline T& operator()(int i, int j, int k) { return mData[index(i, j, k)]; }
125  inline T operator()(int i, int j, int k) const { return mData[index(i, j, k)]; }
127  inline T& operator()(IndexInt idx) { DEBUG_ONLY(checkIndex(idx)); return mData[idx]; }
129  inline T operator()(IndexInt idx) const { DEBUG_ONLY(checkIndex(idx)); return mData[idx]; }
131  inline T& operator()(const Vec3i& pos) { return mData[index(pos)]; }
133  inline T operator()(const Vec3i& pos) const { return mData[index(pos)]; }
135  inline T& operator[](IndexInt idx) { DEBUG_ONLY(checkIndex(idx)); return mData[idx]; }
137  inline const T operator[](IndexInt idx) const { DEBUG_ONLY(checkIndex(idx)); return mData[idx]; }
138 
139  // interpolated access
140  inline T getInterpolated(const Vec3& pos) const { return interpol<T>(mData, mSize, mStrideZ, pos); }
141  inline void setInterpolated(const Vec3& pos, const T& val, Grid<Real>& sumBuffer) const { setInterpol<T>(mData, mSize, mStrideZ, pos, val, &sumBuffer[0]); }
142  // higher order interpolation (1=linear, 2=cubic)
143  inline T getInterpolatedHi(const Vec3& pos, int order) const {
144  switch(order) {
145  case 1: return interpol <T>(mData, mSize, mStrideZ, pos);
146  case 2: return interpolCubic<T>(mData, mSize, mStrideZ, pos);
147  default: assertMsg(false, "Unknown interpolation order "<<order); }
148  }
149 
150  // assignment / copy
151 
153  //Grid<T>& operator=(const Grid<T>& a);
155  PYTHON() Grid<T>& copyFrom(const Grid<T>& a, bool copyType=true ); // old: { *this = a; }
156 
157  // helper functions to work with grids in scene files
158 
160  PYTHON() void add(const Grid<T>& a);
161  PYTHON() void sub(const Grid<T>& a);
163  PYTHON() void setConst(T s);
165  PYTHON() void addConst(T s);
167  PYTHON() void addScaled(const Grid<T>& a, const T& factor);
169  PYTHON() void mult( const Grid<T>& a);
171  PYTHON() void multConst(T s);
173  PYTHON() void clamp(Real min, Real max);
175  PYTHON() void stomp(const T& threshold);
176 
177  // common compound operators
179  PYTHON() Real getMaxAbs();
181  PYTHON() Real getMax();
183  PYTHON() Real getMin();
185  PYTHON() Real getL1(int bnd=0);
187  PYTHON() Real getL2(int bnd=0);
188 
190  PYTHON() void setBound(T value, int boundaryWidth=1);
192  PYTHON() void setBoundNeumann(int boundaryWidth=1);
193 
195  PYTHON() std::string getDataPointer();
196 
198  PYTHON() void printGrid(int zSlice=-1, bool printIndex=false);
199 
200  // c++ only operators
201  template<class S> Grid<T>& operator+=(const Grid<S>& a);
202  template<class S> Grid<T>& operator+=(const S& a);
203  template<class S> Grid<T>& operator-=(const Grid<S>& a);
204  template<class S> Grid<T>& operator-=(const S& a);
205  template<class S> Grid<T>& operator*=(const Grid<S>& a);
206  template<class S> Grid<T>& operator*=(const S& a);
207  template<class S> Grid<T>& operator/=(const Grid<S>& a);
208  template<class S> Grid<T>& operator/=(const S& a);
209  Grid<T>& safeDivide(const Grid<T>& a);
210 
212  void swap(Grid<T>& other);
213 
215  inline T& operator()(int i, int j, int k, int unused) { return mData[index(i, j, k)]; }
216  inline T operator() (int i, int j, int k, int unused) const { return mData[index(i, j, k)]; }
217 
218 protected:
219  T* mData;
220 };
221 
222 // Python doesn't know about templates: explicit aliases needed
223 
224 
225 
226 
229 PYTHON class MACGrid : public Grid<Vec3> {
230 public:
231  PYTHON() MACGrid(FluidSolver* parent, bool show=true) : Grid<Vec3>(parent, show) {
232  mType = (GridType)(TypeMAC | TypeVec3); }
233 
234  // specialized functions for interpolating MAC information
235  inline Vec3 getCentered(int i, int j, int k) const;
236  inline Vec3 getCentered(const Vec3i& pos) const { return getCentered(pos.x, pos.y, pos.z); }
237  inline Vec3 getAtMACX(int i, int j, int k) const;
238  inline Vec3 getAtMACY(int i, int j, int k) const;
239  inline Vec3 getAtMACZ(int i, int j, int k) const;
240  // interpolation
241  inline Vec3 getInterpolated(const Vec3& pos) const { return interpolMAC(mData, mSize, mStrideZ, pos); }
242  inline void setInterpolated(const Vec3& pos, const Vec3& val, Vec3* tmp) { return setInterpolMAC(mData, mSize, mStrideZ, pos, val, tmp); }
243  inline Vec3 getInterpolatedHi(const Vec3& pos, int order) const {
244  switch(order) {
245  case 1: return interpolMAC (mData, mSize, mStrideZ, pos);
246  case 2: return interpolCubicMAC(mData, mSize, mStrideZ, pos);
247  default: assertMsg(false, "Unknown interpolation order "<<order); }
248  }
249  // specials for mac grid:
250  template<int comp> inline Real getInterpolatedComponent(Vec3 pos) const { return interpolComponent<comp>(mData, mSize, mStrideZ, pos); }
251  template<int comp> inline Real getInterpolatedComponentHi(const Vec3& pos, int order) const {
252  switch(order) {
253  case 1: return interpolComponent<comp>(mData, mSize, mStrideZ, pos);
254  case 2: return interpolCubicMAC(mData, mSize, mStrideZ, pos)[comp]; // warning - not yet optimized
255  default: assertMsg(false, "Unknown interpolation order "<<order); }
256  }
257 
260  PYTHON() void setBoundMAC(Vec3 value, int boundaryWidth, bool normalOnly=false);
261 
262 protected:
263 };
264 
267 PYTHON class FlagGrid : public Grid<int> {
268 public:
269  PYTHON() FlagGrid(FluidSolver* parent, int dim=3, bool show=true) : Grid<int>(parent, show) {
270  mType = (GridType)(TypeFlags | TypeInt); }
271 
273  enum CellType {
274  TypeNone = 0,
275  TypeFluid = 1,
276  TypeObstacle = 2,
277  TypeEmpty = 4,
278  TypeInflow = 8,
279  TypeOutflow = 16,
280  TypeOpen = 32,
281  TypeStick = 64,
282  // internal use only, for fast marching
283  TypeReserved = 256,
284  // 2^10 - 2^14 reserved for moving obstacles
285  };
286 
288  inline int getAt(const Vec3& pos) const { return mData[index((int)pos.x, (int)pos.y, (int)pos.z)]; }
289 
291  inline bool isObstacle(IndexInt idx) const { return get(idx) & TypeObstacle; }
292  inline bool isObstacle(int i, int j, int k) const { return get(i,j,k) & TypeObstacle; }
293  inline bool isObstacle(const Vec3i& pos) const { return get(pos) & TypeObstacle; }
294  inline bool isObstacle(const Vec3& pos) const { return getAt(pos) & TypeObstacle; }
295  inline bool isFluid(IndexInt idx) const { return get(idx) & TypeFluid; }
296  inline bool isFluid(int i, int j, int k) const { return get(i,j,k) & TypeFluid; }
297  inline bool isFluid(const Vec3i& pos) const { return get(pos) & TypeFluid; }
298  inline bool isFluid(const Vec3& pos) const { return getAt(pos) & TypeFluid; }
299  inline bool isInflow(IndexInt idx) const { return get(idx) & TypeInflow; }
300  inline bool isInflow(int i, int j, int k) const { return get(i,j,k) & TypeInflow; }
301  inline bool isInflow(const Vec3i& pos) const { return get(pos) & TypeInflow; }
302  inline bool isInflow(const Vec3& pos) const { return getAt(pos) & TypeInflow; }
303  inline bool isEmpty(IndexInt idx) const { return get(idx) & TypeEmpty; }
304  inline bool isEmpty(int i, int j, int k) const { return get(i,j,k) & TypeEmpty; }
305  inline bool isEmpty(const Vec3i& pos) const { return get(pos) & TypeEmpty; }
306  inline bool isEmpty(const Vec3& pos) const { return getAt(pos) & TypeEmpty; }
307  inline bool isOutflow(IndexInt idx) const { return get(idx) & TypeOutflow; }
308  inline bool isOutflow(int i, int j, int k) const { return get(i, j, k) & TypeOutflow; }
309  inline bool isOutflow(const Vec3i& pos) const { return get(pos) & TypeOutflow; }
310  inline bool isOutflow(const Vec3& pos) const { return getAt(pos) & TypeOutflow; }
311  inline bool isOpen(IndexInt idx) const { return get(idx) & TypeOpen; }
312  inline bool isOpen(int i, int j, int k) const { return get(i, j, k) & TypeOpen; }
313  inline bool isOpen(const Vec3i& pos) const { return get(pos) & TypeOpen; }
314  inline bool isOpen(const Vec3& pos) const { return getAt(pos) & TypeOpen; }
315  inline bool isStick(IndexInt idx) const { return get(idx) & TypeStick; }
316  inline bool isStick(int i, int j, int k) const { return get(i,j,k) & TypeStick; }
317  inline bool isStick(const Vec3i& pos) const { return get(pos) & TypeStick; }
318  inline bool isStick(const Vec3& pos) const { return getAt(pos) & TypeStick; }
319 
320 
321  void InitMinXWall(const int &boundaryWidth, Grid<Real>& phiWalls);
322  void InitMaxXWall(const int &boundaryWidth, Grid<Real>& phiWalls);
323  void InitMinYWall(const int &boundaryWidth, Grid<Real>& phiWalls);
324  void InitMaxYWall(const int &boundaryWidth, Grid<Real>& phiWalls);
325  void InitMinZWall(const int &boundaryWidth, Grid<Real>& phiWalls);
326  void InitMaxZWall(const int &boundaryWidth, Grid<Real>& phiWalls);
327  // Python callables
328  PYTHON() void initDomain( const int &boundaryWidth = 0
329  , const std::string &wall = "xXyYzZ"
330  , const std::string &open = " "
331  , const std::string &inflow = " "
332  , const std::string &outflow = " "
333  , Grid<Real>* phiWalls = 0x00 );
334 
335  void initBoundaries( const int &boundaryWidth, const int *types );
336 
338  PYTHON() void updateFromLevelset(LevelsetGrid& levelset);
340  PYTHON() void fillGrid(int type=TypeFluid);
341 
345  PYTHON() int countCells(int flag, int bnd=0, Grid<Real>* mask=NULL);
346 };
347 
350  return Vec3( Real(s1[0])/s2[0], Real(s1[1])/s2[1], Real(s1[2])/s2[2] );
351 }
352 
353 // prototypes for grid plugins
354 void copyMacToVec3(MACGrid &source, Grid<Vec3>& target);
355 void convertMacToVec3(MACGrid &source, Grid<Vec3>& target);
356 void resampleVec3ToMac(Grid<Vec3>& source, MACGrid &target);
357 void resampleMacToVec3 (MACGrid &source, Grid<Vec3>& target );
358 
359 void getComponent(const Grid<Vec3>& source, Grid<Real>& target, int component);
360 void setComponent(const Grid<Real>& source, Grid<Vec3>& target, int component);
361 
362 
363 
364 
365 //******************************************************************************
366 // Implementation of inline functions
367 
368 inline void GridBase::checkIndex(int i, int j, int k) const {
369  //if (i<0 || j<0 || i>=mSize.x || j>=mSize.y || (is3D() && (k<0|| k>= mSize.z))) {
370  if (i<0 || j<0 || i>=mSize.x || j>=mSize.y || k<0|| k>= mSize.z ) {
371  std::ostringstream s;
372  s << "Grid " << mName << " dim " << mSize << " : index " << i << "," << j << "," << k << " out of bound ";
373  errMsg(s.str());
374  }
375 }
376 
377 inline void GridBase::checkIndex(IndexInt idx) const {
378  if (idx<0 || idx >= mSize.x * mSize.y * mSize.z) {
379  std::ostringstream s;
380  s << "Grid " << mName << " dim " << mSize << " : index " << idx << " out of bound ";
381  errMsg(s.str());
382  }
383 }
384 
385 bool GridBase::isInBounds(const Vec3i& p) const {
386  return (p.x >= 0 && p.y >= 0 && p.z >= 0 && p.x < mSize.x && p.y < mSize.y && p.z < mSize.z);
387 }
388 
389 bool GridBase::isInBounds(const Vec3i& p, int bnd) const {
390  bool ret = (p.x >= bnd && p.y >= bnd && p.x < mSize.x-bnd && p.y < mSize.y-bnd);
391  if(this->is3D()) {
392  ret &= (p.z >= bnd && p.z < mSize.z-bnd);
393  } else {
394  ret &= (p.z == 0);
395  }
396  return ret;
397 }
399 bool GridBase::isInBounds(IndexInt idx) const {
400  if (idx<0 || idx >= mSize.x * mSize.y * mSize.z) {
401  return false;
402  }
403  return true;
404 }
405 
406 inline Vec3 MACGrid::getCentered(int i, int j, int k) const {
407  DEBUG_ONLY(checkIndex(i+1,j+1,k));
408  const IndexInt idx = index(i,j,k);
409  Vec3 v = Vec3(0.5* (mData[idx].x + mData[idx+1].x),
410  0.5* (mData[idx].y + mData[idx+mSize.x].y),
411  0.);
412  if( this->is3D() ) {
413  DEBUG_ONLY(checkIndex(idx+mStrideZ));
414  v[2] = 0.5* (mData[idx].z + mData[idx+mStrideZ].z);
415  }
416  return v;
417 }
418 
419 inline Vec3 MACGrid::getAtMACX(int i, int j, int k) const {
420  DEBUG_ONLY(checkIndex(i-1,j+1,k));
421  const IndexInt idx = index(i,j,k);
422  Vec3 v = Vec3( (mData[idx].x),
423  0.25* (mData[idx].y + mData[idx-1].y + mData[idx+mSize.x].y + mData[idx+mSize.x-1].y),
424  0.);
425  if( this->is3D() ) {
426  DEBUG_ONLY(checkIndex(idx+mStrideZ-1));
427  v[2] = 0.25* (mData[idx].z + mData[idx-1].z + mData[idx+mStrideZ].z + mData[idx+mStrideZ-1].z);
428  }
429  return v;
430 }
431 
432 inline Vec3 MACGrid::getAtMACY(int i, int j, int k) const {
433  DEBUG_ONLY(checkIndex(i+1,j-1,k));
434  const IndexInt idx = index(i,j,k);
435  Vec3 v = Vec3(0.25* (mData[idx].x + mData[idx-mSize.x].x + mData[idx+1].x + mData[idx+1-mSize.x].x),
436  (mData[idx].y), 0. );
437  if( this->is3D() ) {
438  DEBUG_ONLY(checkIndex(idx+mStrideZ-mSize.x));
439  v[2] = 0.25* (mData[idx].z + mData[idx-mSize.x].z + mData[idx+mStrideZ].z + mData[idx+mStrideZ-mSize.x].z);
440  }
441  return v;
442 }
443 
444 inline Vec3 MACGrid::getAtMACZ(int i, int j, int k) const {
445  const IndexInt idx = index(i,j,k);
446  DEBUG_ONLY(checkIndex(idx-mStrideZ));
447  DEBUG_ONLY(checkIndex(idx+mSize.x-mStrideZ));
448  Vec3 v = Vec3(0.25* (mData[idx].x + mData[idx-mStrideZ].x + mData[idx+1].x + mData[idx+1-mStrideZ].x),
449  0.25* (mData[idx].y + mData[idx-mStrideZ].y + mData[idx+mSize.x].y + mData[idx+mSize.x-mStrideZ].y),
450  (mData[idx].z) );
451  return v;
452 }
453 
455 void gridAdd (Grid<T>& me, const Grid<S>& other) {}
456 
458 void gridSub (Grid<T>& me, const Grid<S>& other) {}
459 
461 void gridMult (Grid<T>& me, const Grid<S>& other) {}
462 
464 void gridDiv (Grid<T>& me, const Grid<S>& other) {}
465 
467 void gridAddScalar (Grid<T>& me, const S& other) {}
468 
470 void gridMultScalar(Grid<T>& me, const S& other) {}
471 
473 void gridScaledAdd (Grid<T>& me, const Grid<T>& other, const S& factor) {}
474 
475 
477 void gridSetConst(Grid<T>& grid, T value) {}
478 
479 
480 template<class T> template<class S> Grid<T>& Grid<T>::operator+= (const Grid<S>& a) {
481  gridAdd<T,S> (*this, a);
482  return *this;
483 }
484 template<class T> template<class S> Grid<T>& Grid<T>::operator+= (const S& a) {
485  gridAddScalar<T,S> (*this, a);
486  return *this;
487 }
488 template<class T> template<class S> Grid<T>& Grid<T>::operator-= (const Grid<S>& a) {
489  gridSub<T,S> (*this, a);
490  return *this;
491 }
492 template<class T> template<class S> Grid<T>& Grid<T>::operator-= (const S& a) {
493  gridAddScalar<T,S> (*this, -a);
494  return *this;
495 }
496 template<class T> template<class S> Grid<T>& Grid<T>::operator*= (const Grid<S>& a) {
497  gridMult<T,S> (*this, a);
498  return *this;
499 }
500 template<class T> template<class S> Grid<T>& Grid<T>::operator*= (const S& a) {
501  gridMultScalar<T,S> (*this, a);
502  return *this;
503 }
504 template<class T> template<class S> Grid<T>& Grid<T>::operator/= (const Grid<S>& a) {
505  gridDiv<T,S> (*this, a);
506  return *this;
507 }
508 template<class T> template<class S> Grid<T>& Grid<T>::operator/= (const S& a) {
509  S rez((S)1.0 / a);
510  gridMultScalar<T,S> (*this, rez);
511  return *this;
512 }
513 
514 //******************************************************************************
515 // Other helper functions
516 
517 // compute gradient of a scalar grid
518 inline Vec3 getGradient(const Grid<Real>& data, int i, int j, int k) {
519  Vec3 v;
520 
521  if (i > data.getSizeX()-2) i= data.getSizeX()-2;
522  if (j > data.getSizeY()-2) j= data.getSizeY()-2;
523  if (i < 1) i = 1;
524  if (j < 1) j = 1;
525  v = Vec3( data(i+1,j ,k ) - data(i-1,j ,k ) ,
526  data(i ,j+1,k ) - data(i ,j-1,k ) , 0. );
527 
528  if(data.is3D()) {
529  if (k > data.getSizeZ()-2) k= data.getSizeZ()-2;
530  if (k < 1) k = 1;
531  v[2]= data(i ,j ,k+1) - data(i ,j ,k-1);
532  }
533 
534  return v;
535 }
536 
537 // interpolate grid from one size to another size
539 void knInterpolateGridTempl(Grid<S>& target, Grid<S>& source, const Vec3& sourceFactor , Vec3 offset, int orderSpace=1 ) {}
540 
541 // template glue code - choose interpolation based on template arguments
542 template<class GRID>
543 void interpolGridTempl( GRID& target, GRID& source ) {
544  errMsg("interpolGridTempl - Only valid for specific instantiations");
545 }
546 
547 
548 } //namespace
549 #endif
550 
551 
bool isInBounds(const Vec3 &p, int bnd=0) const
Check if index is within given boundaries.
Definition: grid.h:62
Definition: commonkernels.h:22
IndexInt getStrideX() const
Get Stride in X dimension.
Definition: grid.h:45
IndexInt getStrideZ() const
Get Stride in Z dimension.
Definition: grid.h:49
void checkIndex(int i, int j, int k) const
Check if indices are within bounds, otherwise error (should only be called when debugging) ...
Definition: grid.h:368
IndexInt getStrideY() const
Get Stride in Y dimension.
Definition: grid.h:47
bool is4D() const
grid4d compatibility functions
Definition: grid.h:77
int getSizeZ() const
Get the grids Z dimension.
Definition: grid.h:40
Definition: grid.h:267
T operator()(int i, int j, int k) const
access data
Definition: grid.h:125
IndexInt index(const Vec3i &pos) const
Get index into the data.
Definition: grid.h:74
Definition: fileio.h:25
T operator()(IndexInt idx) const
access data
Definition: grid.h:129
Definition: grid.h:229
const T operator[](IndexInt idx) const
access data
Definition: grid.h:137
Basic inlined vector class.
Definition: vectorbase.h:71
T & operator()(IndexInt idx)
access data
Definition: grid.h:127
PYTHON void resampleVec3ToMac(Grid< Vec3 > &source, MACGrid &target)
Definition: grid.cpp:328
bool isObstacle(IndexInt idx) const
check for different flag types
Definition: grid.h:291
int getSizeY() const
Get the grids Y dimension.
Definition: grid.h:38
T & operator()(int i, int j, int k, int unused)
grid4d compatibility functions
Definition: grid.h:215
Vec3i getSize() const
Get the grids dimensions.
Definition: grid.h:42
T & operator()(int i, int j, int k)
access data
Definition: grid.h:123
IndexInt index(int i, int j, int k) const
Get index into the data.
Definition: grid.h:72
int getSizeX() const
Get the grids X dimension.
Definition: grid.h:36
Vec3 calcGridSizeFactor(Vec3i s1, Vec3i s2)
helper to compute grid conversion factor between local coordinates of two grids
Definition: grid.h:349
T & operator[](IndexInt idx)
access data
Definition: grid.h:135
Definition: levelset.h:25
GridType getType() const
Get the type of grid.
Definition: grid.h:67
T & operator()(const Vec3i &pos)
access data
Definition: grid.h:131
T operator()(const Vec3i &pos) const
access data
Definition: grid.h:133
bool isInBounds(const Vec3i &p, int bnd) const
Check if index is within given boundaries.
Definition: grid.h:389
int getAt(const Vec3 &pos) const
access for particles
Definition: grid.h:288
bool is3D() const
Check dimensionality.
Definition: grid.h:69
Definition: grid.h:29
PYTHON void resampleMacToVec3(MACGrid &source, Grid< Vec3 > &target)
Definition: grid.cpp:332
Definition: fluidsolver.h:28