Welcome to mirror list, hosted at ThFree Co, Russian Federation.

git.blender.org/blender.git - Unnamed repository; edit this file 'description' to name the repository.
summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
Diffstat (limited to 'extern/mantaflow/preprocessed/grid.h')
-rw-r--r--extern/mantaflow/preprocessed/grid.h2260
1 files changed, 2260 insertions, 0 deletions
diff --git a/extern/mantaflow/preprocessed/grid.h b/extern/mantaflow/preprocessed/grid.h
new file mode 100644
index 00000000000..bd4e0f99f85
--- /dev/null
+++ b/extern/mantaflow/preprocessed/grid.h
@@ -0,0 +1,2260 @@
+
+
+// DO NOT EDIT !
+// This file is generated using the MantaFlow preprocessor (prep generate).
+
+/******************************************************************************
+ *
+ * MantaFlow fluid solver framework
+ * Copyright 2011 Tobias Pfaff, Nils Thuerey
+ *
+ * This program is free software, distributed under the terms of the
+ * Apache License, Version 2.0
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Grid representation
+ *
+ ******************************************************************************/
+
+#ifndef _GRID_H
+#define _GRID_H
+
+#include "manta.h"
+#include "vectorbase.h"
+#include "interpol.h"
+#include "interpolHigh.h"
+#include "kernel.h"
+
+namespace Manta {
+class LevelsetGrid;
+
+//! Base class for all grids
+class GridBase : public PbClass {
+ public:
+ enum GridType {
+ TypeNone = 0,
+ TypeReal = 1,
+ TypeInt = 2,
+ TypeVec3 = 4,
+ TypeMAC = 8,
+ TypeLevelset = 16,
+ TypeFlags = 32
+ };
+
+ GridBase(FluidSolver *parent);
+ static int _W_0(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ PbClass *obj = Pb::objFromPy(_self);
+ if (obj)
+ delete obj;
+ try {
+ PbArgs _args(_linargs, _kwds);
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(0, "GridBase::GridBase", !noTiming);
+ {
+ ArgLocker _lock;
+ FluidSolver *parent = _args.getPtr<FluidSolver>("parent", 0, &_lock);
+ obj = new GridBase(parent);
+ obj->registerObject(_self, &_args);
+ _args.check();
+ }
+ pbFinalizePlugin(obj->getParent(), "GridBase::GridBase", !noTiming);
+ return 0;
+ }
+ catch (std::exception &e) {
+ pbSetError("GridBase::GridBase", e.what());
+ return -1;
+ }
+ }
+
+ //! Get the grids X dimension
+ inline int getSizeX() const
+ {
+ return mSize.x;
+ }
+ static PyObject *_W_1(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ GridBase *pbo = dynamic_cast<GridBase *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "GridBase::getSizeX", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ pbo->_args.copy(_args);
+ _retval = toPy(pbo->getSizeX());
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "GridBase::getSizeX", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("GridBase::getSizeX", e.what());
+ return 0;
+ }
+ }
+
+ //! Get the grids Y dimension
+ inline int getSizeY() const
+ {
+ return mSize.y;
+ }
+ static PyObject *_W_2(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ GridBase *pbo = dynamic_cast<GridBase *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "GridBase::getSizeY", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ pbo->_args.copy(_args);
+ _retval = toPy(pbo->getSizeY());
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "GridBase::getSizeY", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("GridBase::getSizeY", e.what());
+ return 0;
+ }
+ }
+
+ //! Get the grids Z dimension
+ inline int getSizeZ() const
+ {
+ return mSize.z;
+ }
+ static PyObject *_W_3(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ GridBase *pbo = dynamic_cast<GridBase *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "GridBase::getSizeZ", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ pbo->_args.copy(_args);
+ _retval = toPy(pbo->getSizeZ());
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "GridBase::getSizeZ", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("GridBase::getSizeZ", e.what());
+ return 0;
+ }
+ }
+
+ //! Get the grids dimensions
+ inline Vec3i getSize() const
+ {
+ return mSize;
+ }
+ static PyObject *_W_4(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ GridBase *pbo = dynamic_cast<GridBase *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "GridBase::getSize", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ pbo->_args.copy(_args);
+ _retval = toPy(pbo->getSize());
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "GridBase::getSize", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("GridBase::getSize", e.what());
+ return 0;
+ }
+ }
+
+ //! Get Stride in X dimension
+ inline IndexInt getStrideX() const
+ {
+ return 1;
+ }
+ //! Get Stride in Y dimension
+ inline IndexInt getStrideY() const
+ {
+ return mSize.x;
+ }
+ //! Get Stride in Z dimension
+ inline IndexInt getStrideZ() const
+ {
+ return mStrideZ;
+ }
+
+ inline Real getDx() const
+ {
+ return mDx;
+ }
+
+ //! Check if indices are within bounds, otherwise error (should only be called when debugging)
+ inline void checkIndex(int i, int j, int k) const;
+ //! Check if indices are within bounds, otherwise error (should only be called when debugging)
+ inline void checkIndex(IndexInt idx) const;
+ //! Check if index is within given boundaries
+ inline bool isInBounds(const Vec3i &p, int bnd) const;
+ //! Check if index is within given boundaries
+ inline bool isInBounds(const Vec3i &p) const;
+ //! Check if index is within given boundaries
+ inline bool isInBounds(const Vec3 &p, int bnd = 0) const
+ {
+ return isInBounds(toVec3i(p), bnd);
+ }
+ //! Check if linear index is in the range of the array
+ inline bool isInBounds(IndexInt idx) const;
+
+ //! Get the type of grid
+ inline GridType getType() const
+ {
+ return mType;
+ }
+ //! Check dimensionality
+ inline bool is3D() const
+ {
+ return m3D;
+ }
+ static PyObject *_W_5(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ GridBase *pbo = dynamic_cast<GridBase *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "GridBase::is3D", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ pbo->_args.copy(_args);
+ _retval = toPy(pbo->is3D());
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "GridBase::is3D", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("GridBase::is3D", e.what());
+ return 0;
+ }
+ }
+
+ //! Get index into the data
+ 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;
+ }
+ //! Get index into the data
+ 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;
+ }
+
+ //! grid4d compatibility functions
+ inline bool is4D() const
+ {
+ return false;
+ }
+ static PyObject *_W_6(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ GridBase *pbo = dynamic_cast<GridBase *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "GridBase::is4D", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ pbo->_args.copy(_args);
+ _retval = toPy(pbo->is4D());
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "GridBase::is4D", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("GridBase::is4D", e.what());
+ return 0;
+ }
+ }
+
+ inline int getSizeT() const
+ {
+ return 1;
+ }
+ static PyObject *_W_7(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ GridBase *pbo = dynamic_cast<GridBase *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "GridBase::getSizeT", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ pbo->_args.copy(_args);
+ _retval = toPy(pbo->getSizeT());
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "GridBase::getSizeT", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("GridBase::getSizeT", e.what());
+ return 0;
+ }
+ }
+
+ inline int getStrideT() const
+ {
+ return 0;
+ }
+ static PyObject *_W_8(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ GridBase *pbo = dynamic_cast<GridBase *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "GridBase::getStrideT", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ pbo->_args.copy(_args);
+ _retval = toPy(pbo->getStrideT());
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "GridBase::getStrideT", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("GridBase::getStrideT", e.what());
+ return 0;
+ }
+ }
+
+ inline int index(int i, int j, int k, int unused) const
+ {
+ return index(i, j, k);
+ }
+ 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);
+ }
+
+ protected:
+ GridType mType;
+ Vec3i mSize;
+ Real mDx;
+ bool m3D; // precomputed Z shift: to ensure 2D compatibility, always use this instead of sx*sy !
+ IndexInt mStrideZ;
+ public:
+ PbArgs _args;
+}
+#define _C_GridBase
+;
+
+//! Grid class
+
+template<class T> class Grid : public GridBase {
+ public:
+ //! init new grid, values are set to zero
+ Grid(FluidSolver *parent, bool show = true);
+ static int _W_9(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ PbClass *obj = Pb::objFromPy(_self);
+ if (obj)
+ delete obj;
+ try {
+ PbArgs _args(_linargs, _kwds);
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(0, "Grid::Grid", !noTiming);
+ {
+ ArgLocker _lock;
+ FluidSolver *parent = _args.getPtr<FluidSolver>("parent", 0, &_lock);
+ bool show = _args.getOpt<bool>("show", 1, true, &_lock);
+ obj = new Grid(parent, show);
+ obj->registerObject(_self, &_args);
+ _args.check();
+ }
+ pbFinalizePlugin(obj->getParent(), "Grid::Grid", !noTiming);
+ return 0;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::Grid", e.what());
+ return -1;
+ }
+ }
+
+ //! init new grid with an existing array
+ Grid(FluidSolver *parent, T *data, bool show = true);
+ //! create new & copy content from another grid
+ Grid(const Grid<T> &a);
+ //! return memory to solver
+ virtual ~Grid();
+
+ typedef T BASETYPE;
+ typedef GridBase BASETYPE_GRID;
+
+ void save(std::string name);
+ static PyObject *_W_10(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::save", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ std::string name = _args.get<std::string>("name", 0, &_lock);
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->save(name);
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::save", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::save", e.what());
+ return 0;
+ }
+ }
+
+ void load(std::string name);
+ static PyObject *_W_11(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::load", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ std::string name = _args.get<std::string>("name", 0, &_lock);
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->load(name);
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::load", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::load", e.what());
+ return 0;
+ }
+ }
+
+ //! set all cells to zero
+ void clear();
+ static PyObject *_W_12(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::clear", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->clear();
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::clear", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::clear", e.what());
+ return 0;
+ }
+ }
+
+ //! all kinds of access functions, use grid(), grid[] or grid.get()
+ //! access data
+ inline T get(int i, int j, int k) const
+ {
+ return mData[index(i, j, k)];
+ }
+ //! access data
+ inline T &get(int i, int j, int k)
+ {
+ return mData[index(i, j, k)];
+ }
+ //! access data
+ inline T get(IndexInt idx) const
+ {
+ DEBUG_ONLY(checkIndex(idx));
+ return mData[idx];
+ }
+ //! access data
+ inline T get(const Vec3i &pos) const
+ {
+ return mData[index(pos)];
+ }
+ //! access data
+ inline T &operator()(int i, int j, int k)
+ {
+ return mData[index(i, j, k)];
+ }
+ //! access data
+ inline T operator()(int i, int j, int k) const
+ {
+ return mData[index(i, j, k)];
+ }
+ //! access data
+ inline T &operator()(IndexInt idx)
+ {
+ DEBUG_ONLY(checkIndex(idx));
+ return mData[idx];
+ }
+ //! access data
+ inline T operator()(IndexInt idx) const
+ {
+ DEBUG_ONLY(checkIndex(idx));
+ return mData[idx];
+ }
+ //! access data
+ inline T &operator()(const Vec3i &pos)
+ {
+ return mData[index(pos)];
+ }
+ //! access data
+ inline T operator()(const Vec3i &pos) const
+ {
+ return mData[index(pos)];
+ }
+ //! access data
+ inline T &operator[](IndexInt idx)
+ {
+ DEBUG_ONLY(checkIndex(idx));
+ return mData[idx];
+ }
+ //! access data
+ inline const T operator[](IndexInt idx) const
+ {
+ DEBUG_ONLY(checkIndex(idx));
+ return mData[idx];
+ }
+
+ // interpolated access
+ inline T getInterpolated(const Vec3 &pos) const
+ {
+ return interpol<T>(mData, mSize, mStrideZ, pos);
+ }
+ inline void setInterpolated(const Vec3 &pos, const T &val, Grid<Real> &sumBuffer) const
+ {
+ setInterpol<T>(mData, mSize, mStrideZ, pos, val, &sumBuffer[0]);
+ }
+ // higher order interpolation (1=linear, 2=cubic)
+ inline T getInterpolatedHi(const Vec3 &pos, int order) const
+ {
+ switch (order) {
+ case 1:
+ return interpol<T>(mData, mSize, mStrideZ, pos);
+ case 2:
+ return interpolCubic<T>(mData, mSize, mStrideZ, pos);
+ default:
+ assertMsg(false, "Unknown interpolation order " << order);
+ }
+ return T(0.); // should never be reached, just to prevent compiler warnings
+ }
+
+ // assignment / copy
+
+ //! warning - do not use "=" for grids in python, this copies the reference! not the grid
+ //! content...
+ // Grid<T>& operator=(const Grid<T>& a);
+ //! copy content from other grid (use this one instead of operator= !)
+ Grid<T> &copyFrom(const Grid<T> &a, bool copyType = true);
+ static PyObject *_W_13(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::copyFrom", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ const Grid<T> &a = *_args.getPtr<Grid<T>>("a", 0, &_lock);
+ bool copyType = _args.getOpt<bool>("copyType", 1, true, &_lock);
+ pbo->_args.copy(_args);
+ _retval = toPy(pbo->copyFrom(a, copyType));
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::copyFrom", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::copyFrom", e.what());
+ return 0;
+ }
+ }
+ // old: { *this = a; }
+
+ // helper functions to work with grids in scene files
+
+ //! get grid type
+ int getGridType();
+ static PyObject *_W_14(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::getGridType", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ pbo->_args.copy(_args);
+ _retval = toPy(pbo->getGridType());
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::getGridType", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::getGridType", e.what());
+ return 0;
+ }
+ }
+
+ //! add/subtract other grid
+ void add(const Grid<T> &a);
+ static PyObject *_W_15(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::add", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ const Grid<T> &a = *_args.getPtr<Grid<T>>("a", 0, &_lock);
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->add(a);
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::add", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::add", e.what());
+ return 0;
+ }
+ }
+
+ void sub(const Grid<T> &a);
+ static PyObject *_W_16(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::sub", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ const Grid<T> &a = *_args.getPtr<Grid<T>>("a", 0, &_lock);
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->sub(a);
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::sub", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::sub", e.what());
+ return 0;
+ }
+ }
+
+ //! set all cells to constant value
+ void setConst(T s);
+ static PyObject *_W_17(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::setConst", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ T s = _args.get<T>("s", 0, &_lock);
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->setConst(s);
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::setConst", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::setConst", e.what());
+ return 0;
+ }
+ }
+
+ //! add constant to all grid cells
+ void addConst(T s);
+ static PyObject *_W_18(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::addConst", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ T s = _args.get<T>("s", 0, &_lock);
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->addConst(s);
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::addConst", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::addConst", e.what());
+ return 0;
+ }
+ }
+
+ //! add scaled other grid to current one (note, only "Real" factor, "T" type not supported here!)
+ void addScaled(const Grid<T> &a, const T &factor);
+ static PyObject *_W_19(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::addScaled", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ const Grid<T> &a = *_args.getPtr<Grid<T>>("a", 0, &_lock);
+ const T &factor = *_args.getPtr<T>("factor", 1, &_lock);
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->addScaled(a, factor);
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::addScaled", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::addScaled", e.what());
+ return 0;
+ }
+ }
+
+ //! multiply contents of grid
+ void mult(const Grid<T> &a);
+ static PyObject *_W_20(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::mult", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ const Grid<T> &a = *_args.getPtr<Grid<T>>("a", 0, &_lock);
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->mult(a);
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::mult", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::mult", e.what());
+ return 0;
+ }
+ }
+
+ //! multiply each cell by a constant scalar value
+ void multConst(T s);
+ static PyObject *_W_21(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::multConst", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ T s = _args.get<T>("s", 0, &_lock);
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->multConst(s);
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::multConst", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::multConst", e.what());
+ return 0;
+ }
+ }
+
+ //! clamp content to range (for vec3, clamps each component separately)
+ void clamp(Real min, Real max);
+ static PyObject *_W_22(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::clamp", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ Real min = _args.get<Real>("min", 0, &_lock);
+ Real max = _args.get<Real>("max", 1, &_lock);
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->clamp(min, max);
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::clamp", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::clamp", e.what());
+ return 0;
+ }
+ }
+
+ //! reduce small values to zero
+ void stomp(const T &threshold);
+ static PyObject *_W_23(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::stomp", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ const T &threshold = *_args.getPtr<T>("threshold", 0, &_lock);
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->stomp(threshold);
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::stomp", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::stomp", e.what());
+ return 0;
+ }
+ }
+
+ //! permute grid axes, e.g. switch y with z (0,2,1)
+ void permuteAxes(int axis0, int axis1, int axis2);
+ static PyObject *_W_24(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::permuteAxes", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ int axis0 = _args.get<int>("axis0", 0, &_lock);
+ int axis1 = _args.get<int>("axis1", 1, &_lock);
+ int axis2 = _args.get<int>("axis2", 2, &_lock);
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->permuteAxes(axis0, axis1, axis2);
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::permuteAxes", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::permuteAxes", e.what());
+ return 0;
+ }
+ }
+
+ //! permute grid axes, e.g. switch y with z (0,2,1)
+ void permuteAxesCopyToGrid(int axis0, int axis1, int axis2, Grid<T> &out);
+ static PyObject *_W_25(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::permuteAxesCopyToGrid", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ int axis0 = _args.get<int>("axis0", 0, &_lock);
+ int axis1 = _args.get<int>("axis1", 1, &_lock);
+ int axis2 = _args.get<int>("axis2", 2, &_lock);
+ Grid<T> &out = *_args.getPtr<Grid<T>>("out", 3, &_lock);
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->permuteAxesCopyToGrid(axis0, axis1, axis2, out);
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::permuteAxesCopyToGrid", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::permuteAxesCopyToGrid", e.what());
+ return 0;
+ }
+ }
+
+ // common compound operators
+ //! get absolute max value in grid
+ Real getMaxAbs() const;
+ static PyObject *_W_26(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::getMaxAbs", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ pbo->_args.copy(_args);
+ _retval = toPy(pbo->getMaxAbs());
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::getMaxAbs", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::getMaxAbs", e.what());
+ return 0;
+ }
+ }
+
+ //! get max value in grid
+ Real getMax() const;
+ static PyObject *_W_27(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::getMax", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ pbo->_args.copy(_args);
+ _retval = toPy(pbo->getMax());
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::getMax", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::getMax", e.what());
+ return 0;
+ }
+ }
+
+ //! get min value in grid
+ Real getMin() const;
+ static PyObject *_W_28(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::getMin", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ pbo->_args.copy(_args);
+ _retval = toPy(pbo->getMin());
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::getMin", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::getMin", e.what());
+ return 0;
+ }
+ }
+
+ //! calculate L1 norm of grid content
+ Real getL1(int bnd = 0);
+ static PyObject *_W_29(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::getL1", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ int bnd = _args.getOpt<int>("bnd", 0, 0, &_lock);
+ pbo->_args.copy(_args);
+ _retval = toPy(pbo->getL1(bnd));
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::getL1", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::getL1", e.what());
+ return 0;
+ }
+ }
+
+ //! calculate L2 norm of grid content
+ Real getL2(int bnd = 0);
+ static PyObject *_W_30(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::getL2", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ int bnd = _args.getOpt<int>("bnd", 0, 0, &_lock);
+ pbo->_args.copy(_args);
+ _retval = toPy(pbo->getL2(bnd));
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::getL2", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::getL2", e.what());
+ return 0;
+ }
+ }
+
+ //! set all boundary cells to constant value (Dirichlet)
+ void setBound(T value, int boundaryWidth = 1);
+ static PyObject *_W_31(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::setBound", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ T value = _args.get<T>("value", 0, &_lock);
+ int boundaryWidth = _args.getOpt<int>("boundaryWidth", 1, 1, &_lock);
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->setBound(value, boundaryWidth);
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::setBound", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::setBound", e.what());
+ return 0;
+ }
+ }
+
+ //! set all boundary cells to last inner value (Neumann)
+ void setBoundNeumann(int boundaryWidth = 1);
+ static PyObject *_W_32(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::setBoundNeumann", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ int boundaryWidth = _args.getOpt<int>("boundaryWidth", 0, 1, &_lock);
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->setBoundNeumann(boundaryWidth);
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::setBoundNeumann", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::setBoundNeumann", e.what());
+ return 0;
+ }
+ }
+
+ //! get data pointer of grid
+ std::string getDataPointer();
+ static PyObject *_W_33(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::getDataPointer", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ pbo->_args.copy(_args);
+ _retval = toPy(pbo->getDataPointer());
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::getDataPointer", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::getDataPointer", e.what());
+ return 0;
+ }
+ }
+
+ //! debugging helper, print grid from python. skip boundary of width bnd
+ void printGrid(int zSlice = -1, bool printIndex = false, int bnd = 1);
+ static PyObject *_W_34(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ Grid *pbo = dynamic_cast<Grid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "Grid::printGrid", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ int zSlice = _args.getOpt<int>("zSlice", 0, -1, &_lock);
+ bool printIndex = _args.getOpt<bool>("printIndex", 1, false, &_lock);
+ int bnd = _args.getOpt<int>("bnd", 2, 1, &_lock);
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->printGrid(zSlice, printIndex, bnd);
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "Grid::printGrid", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("Grid::printGrid", e.what());
+ return 0;
+ }
+ }
+
+ // c++ only operators
+ template<class S> Grid<T> &operator+=(const Grid<S> &a);
+ template<class S> Grid<T> &operator+=(const S &a);
+ template<class S> Grid<T> &operator-=(const Grid<S> &a);
+ template<class S> Grid<T> &operator-=(const S &a);
+ template<class S> Grid<T> &operator*=(const Grid<S> &a);
+ template<class S> Grid<T> &operator*=(const S &a);
+ template<class S> Grid<T> &operator/=(const Grid<S> &a);
+ template<class S> Grid<T> &operator/=(const S &a);
+ Grid<T> &safeDivide(const Grid<T> &a);
+
+ //! Swap data with another grid (no actual data is moved)
+ void swap(Grid<T> &other);
+
+ //! grid4d compatibility functions
+ inline T &operator()(int i, int j, int k, int unused)
+ {
+ return mData[index(i, j, k)];
+ }
+ inline T operator()(int i, int j, int k, int unused) const
+ {
+ return mData[index(i, j, k)];
+ }
+
+ protected:
+ T *mData;
+ bool externalData; // True if mData is managed outside of the Fluidsolver
+ public:
+ PbArgs _args;
+}
+#define _C_Grid
+;
+
+// Python doesn't know about templates: explicit aliases needed
+
+//! Special function for staggered grids
+class MACGrid : public Grid<Vec3> {
+ public:
+ MACGrid(FluidSolver *parent, bool show = true) : Grid<Vec3>(parent, show)
+ {
+ mType = (GridType)(TypeMAC | TypeVec3);
+ }
+ static int _W_35(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ PbClass *obj = Pb::objFromPy(_self);
+ if (obj)
+ delete obj;
+ try {
+ PbArgs _args(_linargs, _kwds);
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(0, "MACGrid::MACGrid", !noTiming);
+ {
+ ArgLocker _lock;
+ FluidSolver *parent = _args.getPtr<FluidSolver>("parent", 0, &_lock);
+ bool show = _args.getOpt<bool>("show", 1, true, &_lock);
+ obj = new MACGrid(parent, show);
+ obj->registerObject(_self, &_args);
+ _args.check();
+ }
+ pbFinalizePlugin(obj->getParent(), "MACGrid::MACGrid", !noTiming);
+ return 0;
+ }
+ catch (std::exception &e) {
+ pbSetError("MACGrid::MACGrid", e.what());
+ return -1;
+ }
+ }
+
+ MACGrid(FluidSolver *parent, Vec3 *data, bool show = true) : Grid<Vec3>(parent, data, show)
+ {
+ mType = (GridType)(TypeMAC | TypeVec3);
+ }
+
+ // specialized functions for interpolating MAC information
+ inline Vec3 getCentered(int i, int j, int k) const;
+ inline Vec3 getCentered(const Vec3i &pos) const
+ {
+ return getCentered(pos.x, pos.y, pos.z);
+ }
+ inline Vec3 getAtMACX(int i, int j, int k) const;
+ inline Vec3 getAtMACY(int i, int j, int k) const;
+ inline Vec3 getAtMACZ(int i, int j, int k) const;
+ // interpolation
+ inline Vec3 getInterpolated(const Vec3 &pos) const
+ {
+ return interpolMAC(mData, mSize, mStrideZ, pos);
+ }
+ inline void setInterpolated(const Vec3 &pos, const Vec3 &val, Vec3 *tmp) const
+ {
+ return setInterpolMAC(mData, mSize, mStrideZ, pos, val, tmp);
+ }
+ inline Vec3 getInterpolatedHi(const Vec3 &pos, int order) const
+ {
+ switch (order) {
+ case 1:
+ return interpolMAC(mData, mSize, mStrideZ, pos);
+ case 2:
+ return interpolCubicMAC(mData, mSize, mStrideZ, pos);
+ default:
+ assertMsg(false, "Unknown interpolation order " << order);
+ }
+ return Vec3(0.); // should never be reached, just to prevent compiler warnings
+ }
+ // specials for mac grid:
+ template<int comp> inline Real getInterpolatedComponent(Vec3 pos) const
+ {
+ return interpolComponent<comp>(mData, mSize, mStrideZ, pos);
+ }
+ template<int comp> inline Real getInterpolatedComponentHi(const Vec3 &pos, int order) const
+ {
+ switch (order) {
+ case 1:
+ return interpolComponent<comp>(mData, mSize, mStrideZ, pos);
+ case 2:
+ return interpolCubicMAC(mData, mSize, mStrideZ, pos)[comp]; // warning - not yet optimized
+ default:
+ assertMsg(false, "Unknown interpolation order " << order);
+ }
+ return 0.; // should never be reached, just to prevent compiler warnings
+ }
+
+ //! set all boundary cells of a MAC grid to certain value (Dirchlet). Respects staggered grid
+ //! locations optionally, only set normal components
+ void setBoundMAC(Vec3 value, int boundaryWidth, bool normalOnly = false);
+ static PyObject *_W_36(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ MACGrid *pbo = dynamic_cast<MACGrid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "MACGrid::setBoundMAC", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ Vec3 value = _args.get<Vec3>("value", 0, &_lock);
+ int boundaryWidth = _args.get<int>("boundaryWidth", 1, &_lock);
+ bool normalOnly = _args.getOpt<bool>("normalOnly", 2, false, &_lock);
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->setBoundMAC(value, boundaryWidth, normalOnly);
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "MACGrid::setBoundMAC", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("MACGrid::setBoundMAC", e.what());
+ return 0;
+ }
+ }
+
+ protected:
+ public:
+ PbArgs _args;
+}
+#define _C_MACGrid
+;
+
+//! Special functions for FlagGrid
+class FlagGrid : public Grid<int> {
+ public:
+ FlagGrid(FluidSolver *parent, int dim = 3, bool show = true) : Grid<int>(parent, show)
+ {
+ mType = (GridType)(TypeFlags | TypeInt);
+ }
+ static int _W_37(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ PbClass *obj = Pb::objFromPy(_self);
+ if (obj)
+ delete obj;
+ try {
+ PbArgs _args(_linargs, _kwds);
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(0, "FlagGrid::FlagGrid", !noTiming);
+ {
+ ArgLocker _lock;
+ FluidSolver *parent = _args.getPtr<FluidSolver>("parent", 0, &_lock);
+ int dim = _args.getOpt<int>("dim", 1, 3, &_lock);
+ bool show = _args.getOpt<bool>("show", 2, true, &_lock);
+ obj = new FlagGrid(parent, dim, show);
+ obj->registerObject(_self, &_args);
+ _args.check();
+ }
+ pbFinalizePlugin(obj->getParent(), "FlagGrid::FlagGrid", !noTiming);
+ return 0;
+ }
+ catch (std::exception &e) {
+ pbSetError("FlagGrid::FlagGrid", e.what());
+ return -1;
+ }
+ }
+
+ FlagGrid(FluidSolver *parent, int *data, int dim = 3, bool show = true)
+ : Grid<int>(parent, data, show)
+ {
+ mType = (GridType)(TypeFlags | TypeInt);
+ }
+
+ //! types of cells, in/outflow can be combined, e.g., TypeFluid|TypeInflow
+ enum CellType {
+ TypeNone = 0,
+ TypeFluid = 1,
+ TypeObstacle = 2,
+ TypeEmpty = 4,
+ TypeInflow = 8,
+ TypeOutflow = 16,
+ TypeOpen = 32,
+ TypeStick = 64,
+ TypeReserved = 256
+ };
+
+ //! access for particles
+ inline int getAt(const Vec3 &pos) const
+ {
+ return mData[index((int)pos.x, (int)pos.y, (int)pos.z)];
+ }
+
+ //! check for different flag types
+ inline bool isObstacle(IndexInt idx) const
+ {
+ return get(idx) & TypeObstacle;
+ }
+ inline bool isObstacle(int i, int j, int k) const
+ {
+ return get(i, j, k) & TypeObstacle;
+ }
+ inline bool isObstacle(const Vec3i &pos) const
+ {
+ return get(pos) & TypeObstacle;
+ }
+ inline bool isObstacle(const Vec3 &pos) const
+ {
+ return getAt(pos) & TypeObstacle;
+ }
+ inline bool isFluid(IndexInt idx) const
+ {
+ return get(idx) & TypeFluid;
+ }
+ inline bool isFluid(int i, int j, int k) const
+ {
+ return get(i, j, k) & TypeFluid;
+ }
+ inline bool isFluid(const Vec3i &pos) const
+ {
+ return get(pos) & TypeFluid;
+ }
+ inline bool isFluid(const Vec3 &pos) const
+ {
+ return getAt(pos) & TypeFluid;
+ }
+ inline bool isInflow(IndexInt idx) const
+ {
+ return get(idx) & TypeInflow;
+ }
+ inline bool isInflow(int i, int j, int k) const
+ {
+ return get(i, j, k) & TypeInflow;
+ }
+ inline bool isInflow(const Vec3i &pos) const
+ {
+ return get(pos) & TypeInflow;
+ }
+ inline bool isInflow(const Vec3 &pos) const
+ {
+ return getAt(pos) & TypeInflow;
+ }
+ inline bool isEmpty(IndexInt idx) const
+ {
+ return get(idx) & TypeEmpty;
+ }
+ inline bool isEmpty(int i, int j, int k) const
+ {
+ return get(i, j, k) & TypeEmpty;
+ }
+ inline bool isEmpty(const Vec3i &pos) const
+ {
+ return get(pos) & TypeEmpty;
+ }
+ inline bool isEmpty(const Vec3 &pos) const
+ {
+ return getAt(pos) & TypeEmpty;
+ }
+ inline bool isOutflow(IndexInt idx) const
+ {
+ return get(idx) & TypeOutflow;
+ }
+ inline bool isOutflow(int i, int j, int k) const
+ {
+ return get(i, j, k) & TypeOutflow;
+ }
+ inline bool isOutflow(const Vec3i &pos) const
+ {
+ return get(pos) & TypeOutflow;
+ }
+ inline bool isOutflow(const Vec3 &pos) const
+ {
+ return getAt(pos) & TypeOutflow;
+ }
+ inline bool isOpen(IndexInt idx) const
+ {
+ return get(idx) & TypeOpen;
+ }
+ inline bool isOpen(int i, int j, int k) const
+ {
+ return get(i, j, k) & TypeOpen;
+ }
+ inline bool isOpen(const Vec3i &pos) const
+ {
+ return get(pos) & TypeOpen;
+ }
+ inline bool isOpen(const Vec3 &pos) const
+ {
+ return getAt(pos) & TypeOpen;
+ }
+ inline bool isStick(IndexInt idx) const
+ {
+ return get(idx) & TypeStick;
+ }
+ inline bool isStick(int i, int j, int k) const
+ {
+ return get(i, j, k) & TypeStick;
+ }
+ inline bool isStick(const Vec3i &pos) const
+ {
+ return get(pos) & TypeStick;
+ }
+ inline bool isStick(const Vec3 &pos) const
+ {
+ return getAt(pos) & TypeStick;
+ }
+
+ void InitMinXWall(const int &boundaryWidth, Grid<Real> &phiWalls);
+ void InitMaxXWall(const int &boundaryWidth, Grid<Real> &phiWalls);
+ void InitMinYWall(const int &boundaryWidth, Grid<Real> &phiWalls);
+ void InitMaxYWall(const int &boundaryWidth, Grid<Real> &phiWalls);
+ void InitMinZWall(const int &boundaryWidth, Grid<Real> &phiWalls);
+ void InitMaxZWall(const int &boundaryWidth, Grid<Real> &phiWalls);
+ // Python callables
+
+ void initDomain(const int &boundaryWidth = 0,
+ const std::string &wall = "xXyYzZ",
+ const std::string &open = " ",
+ const std::string &inflow = " ",
+ const std::string &outflow = " ",
+ Grid<Real> *phiWalls = 0x00);
+ static PyObject *_W_38(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ FlagGrid *pbo = dynamic_cast<FlagGrid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "FlagGrid::initDomain", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ const int &boundaryWidth = _args.getOpt<int>("boundaryWidth", 0, 0, &_lock);
+ const std::string &wall = _args.getOpt<std::string>("wall", 1, "xXyYzZ", &_lock);
+ const std::string &open = _args.getOpt<std::string>("open", 2, " ", &_lock);
+ const std::string &inflow = _args.getOpt<std::string>("inflow", 3, " ", &_lock);
+ const std::string &outflow = _args.getOpt<std::string>("outflow", 4, " ", &_lock);
+ Grid<Real> *phiWalls = _args.getPtrOpt<Grid<Real>>("phiWalls", 5, 0x00, &_lock);
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->initDomain(boundaryWidth, wall, open, inflow, outflow, phiWalls);
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "FlagGrid::initDomain", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("FlagGrid::initDomain", e.what());
+ return 0;
+ }
+ }
+
+ void initBoundaries(const int &boundaryWidth, const int *types);
+
+ //! set fluid flags inside levelset (liquids)
+ void updateFromLevelset(LevelsetGrid &levelset);
+ static PyObject *_W_39(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ FlagGrid *pbo = dynamic_cast<FlagGrid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "FlagGrid::updateFromLevelset", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ LevelsetGrid &levelset = *_args.getPtr<LevelsetGrid>("levelset", 0, &_lock);
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->updateFromLevelset(levelset);
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "FlagGrid::updateFromLevelset", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("FlagGrid::updateFromLevelset", e.what());
+ return 0;
+ }
+ }
+
+ //! set all cells (except obs/in/outflow) to type (fluid by default)
+ void fillGrid(int type = TypeFluid);
+ static PyObject *_W_40(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ FlagGrid *pbo = dynamic_cast<FlagGrid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "FlagGrid::fillGrid", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ int type = _args.getOpt<int>("type", 0, TypeFluid, &_lock);
+ pbo->_args.copy(_args);
+ _retval = getPyNone();
+ pbo->fillGrid(type);
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "FlagGrid::fillGrid", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("FlagGrid::fillGrid", e.what());
+ return 0;
+ }
+ }
+
+ //! count no. of cells matching flags via "AND"
+ //! warning for large grids! only regular int returned (due to python interface)
+ //! optionally creates mask in RealGrid (1 where flag matches, 0 otherwise)
+ int countCells(int flag, int bnd = 0, Grid<Real> *mask = NULL);
+ static PyObject *_W_41(PyObject *_self, PyObject *_linargs, PyObject *_kwds)
+ {
+ try {
+ PbArgs _args(_linargs, _kwds);
+ FlagGrid *pbo = dynamic_cast<FlagGrid *>(Pb::objFromPy(_self));
+ bool noTiming = _args.getOpt<bool>("notiming", -1, 0);
+ pbPreparePlugin(pbo->getParent(), "FlagGrid::countCells", !noTiming);
+ PyObject *_retval = 0;
+ {
+ ArgLocker _lock;
+ int flag = _args.get<int>("flag", 0, &_lock);
+ int bnd = _args.getOpt<int>("bnd", 1, 0, &_lock);
+ Grid<Real> *mask = _args.getPtrOpt<Grid<Real>>("mask", 2, NULL, &_lock);
+ pbo->_args.copy(_args);
+ _retval = toPy(pbo->countCells(flag, bnd, mask));
+ pbo->_args.check();
+ }
+ pbFinalizePlugin(pbo->getParent(), "FlagGrid::countCells", !noTiming);
+ return _retval;
+ }
+ catch (std::exception &e) {
+ pbSetError("FlagGrid::countCells", e.what());
+ return 0;
+ }
+ }
+
+ public:
+ PbArgs _args;
+}
+#define _C_FlagGrid
+;
+
+//! helper to compute grid conversion factor between local coordinates of two grids
+inline Vec3 calcGridSizeFactor(Vec3i s1, Vec3i s2)
+{
+ return Vec3(Real(s1[0]) / s2[0], Real(s1[1]) / s2[1], Real(s1[2]) / s2[2]);
+}
+
+// prototypes for grid plugins
+void copyMacToVec3(MACGrid &source, Grid<Vec3> &target);
+void convertMacToVec3(MACGrid &source, Grid<Vec3> &target);
+void resampleVec3ToMac(Grid<Vec3> &source, MACGrid &target);
+void resampleMacToVec3(MACGrid &source, Grid<Vec3> &target);
+
+void getComponent(const Grid<Vec3> &source, Grid<Real> &target, int component);
+void setComponent(const Grid<Real> &source, Grid<Vec3> &target, int component);
+
+//******************************************************************************
+// Implementation of inline functions
+
+inline void GridBase::checkIndex(int i, int j, int k) const
+{
+ if (i < 0 || j < 0 || k < 0 || i >= mSize.x || j >= mSize.y || k >= mSize.z) {
+ std::ostringstream s;
+ s << "Grid " << mName << " dim " << mSize << " : index " << i << "," << j << "," << k
+ << " out of bound ";
+ errMsg(s.str());
+ }
+}
+
+inline void GridBase::checkIndex(IndexInt idx) const
+{
+ if (idx < 0 || idx >= mSize.x * mSize.y * mSize.z) {
+ std::ostringstream s;
+ s << "Grid " << mName << " dim " << mSize << " : index " << idx << " out of bound ";
+ errMsg(s.str());
+ }
+}
+
+bool GridBase::isInBounds(const Vec3i &p) const
+{
+ return (p.x >= 0 && p.y >= 0 && p.z >= 0 && p.x < mSize.x && p.y < mSize.y && p.z < mSize.z);
+}
+
+bool GridBase::isInBounds(const Vec3i &p, int bnd) const
+{
+ bool ret = (p.x >= bnd && p.y >= bnd && p.x < mSize.x - bnd && p.y < mSize.y - bnd);
+ if (this->is3D()) {
+ ret &= (p.z >= bnd && p.z < mSize.z - bnd);
+ }
+ else {
+ ret &= (p.z == 0);
+ }
+ return ret;
+}
+//! Check if linear index is in the range of the array
+bool GridBase::isInBounds(IndexInt idx) const
+{
+ if (idx < 0 || idx >= mSize.x * mSize.y * mSize.z) {
+ return false;
+ }
+ return true;
+}
+
+inline Vec3 MACGrid::getCentered(int i, int j, int k) const
+{
+ DEBUG_ONLY(checkIndex(i + 1, j + 1, k));
+ const IndexInt idx = index(i, j, k);
+ Vec3 v = Vec3(
+ 0.5 * (mData[idx].x + mData[idx + 1].x), 0.5 * (mData[idx].y + mData[idx + mSize.x].y), 0.);
+ if (this->is3D()) {
+ DEBUG_ONLY(checkIndex(idx + mStrideZ));
+ v[2] = 0.5 * (mData[idx].z + mData[idx + mStrideZ].z);
+ }
+ return v;
+}
+
+inline Vec3 MACGrid::getAtMACX(int i, int j, int k) const
+{
+ DEBUG_ONLY(checkIndex(i - 1, j + 1, k));
+ const IndexInt idx = index(i, j, k);
+ Vec3 v = Vec3((mData[idx].x),
+ 0.25 * (mData[idx].y + mData[idx - 1].y + mData[idx + mSize.x].y +
+ mData[idx + mSize.x - 1].y),
+ 0.);
+ if (this->is3D()) {
+ DEBUG_ONLY(checkIndex(idx + mStrideZ - 1));
+ v[2] = 0.25 * (mData[idx].z + mData[idx - 1].z + mData[idx + mStrideZ].z +
+ mData[idx + mStrideZ - 1].z);
+ }
+ return v;
+}
+
+inline Vec3 MACGrid::getAtMACY(int i, int j, int k) const
+{
+ DEBUG_ONLY(checkIndex(i + 1, j - 1, k));
+ const IndexInt idx = index(i, j, k);
+ Vec3 v = Vec3(0.25 * (mData[idx].x + mData[idx - mSize.x].x + mData[idx + 1].x +
+ mData[idx + 1 - mSize.x].x),
+ (mData[idx].y),
+ 0.);
+ if (this->is3D()) {
+ DEBUG_ONLY(checkIndex(idx + mStrideZ - mSize.x));
+ v[2] = 0.25 * (mData[idx].z + mData[idx - mSize.x].z + mData[idx + mStrideZ].z +
+ mData[idx + mStrideZ - mSize.x].z);
+ }
+ return v;
+}
+
+inline Vec3 MACGrid::getAtMACZ(int i, int j, int k) const
+{
+ const IndexInt idx = index(i, j, k);
+ DEBUG_ONLY(checkIndex(idx - mStrideZ));
+ DEBUG_ONLY(checkIndex(idx + mSize.x - mStrideZ));
+ Vec3 v = Vec3(0.25 * (mData[idx].x + mData[idx - mStrideZ].x + mData[idx + 1].x +
+ mData[idx + 1 - mStrideZ].x),
+ 0.25 * (mData[idx].y + mData[idx - mStrideZ].y + mData[idx + mSize.x].y +
+ mData[idx + mSize.x - mStrideZ].y),
+ (mData[idx].z));
+ return v;
+}
+
+template<class T, class S> struct gridAdd : public KernelBase {
+ gridAdd(Grid<T> &me, const Grid<S> &other) : KernelBase(&me, 0), me(me), other(other)
+ {
+ runMessage();
+ run();
+ }
+ inline void op(IndexInt idx, Grid<T> &me, const Grid<S> &other) const
+ {
+ me[idx] += other[idx];
+ }
+ inline Grid<T> &getArg0()
+ {
+ return me;
+ }
+ typedef Grid<T> type0;
+ inline const Grid<S> &getArg1()
+ {
+ return other;
+ }
+ typedef Grid<S> type1;
+ void runMessage()
+ {
+ debMsg("Executing kernel gridAdd ", 3);
+ debMsg("Kernel range"
+ << " x " << maxX << " y " << maxY << " z " << minZ << " - " << maxZ << " ",
+ 4);
+ };
+ void operator()(const tbb::blocked_range<IndexInt> &__r) const
+ {
+ for (IndexInt idx = __r.begin(); idx != (IndexInt)__r.end(); idx++)
+ op(idx, me, other);
+ }
+ void run()
+ {
+ tbb::parallel_for(tbb::blocked_range<IndexInt>(0, size), *this);
+ }
+ Grid<T> &me;
+ const Grid<S> &other;
+};
+template<class T, class S> struct gridSub : public KernelBase {
+ gridSub(Grid<T> &me, const Grid<S> &other) : KernelBase(&me, 0), me(me), other(other)
+ {
+ runMessage();
+ run();
+ }
+ inline void op(IndexInt idx, Grid<T> &me, const Grid<S> &other) const
+ {
+ me[idx] -= other[idx];
+ }
+ inline Grid<T> &getArg0()
+ {
+ return me;
+ }
+ typedef Grid<T> type0;
+ inline const Grid<S> &getArg1()
+ {
+ return other;
+ }
+ typedef Grid<S> type1;
+ void runMessage()
+ {
+ debMsg("Executing kernel gridSub ", 3);
+ debMsg("Kernel range"
+ << " x " << maxX << " y " << maxY << " z " << minZ << " - " << maxZ << " ",
+ 4);
+ };
+ void operator()(const tbb::blocked_range<IndexInt> &__r) const
+ {
+ for (IndexInt idx = __r.begin(); idx != (IndexInt)__r.end(); idx++)
+ op(idx, me, other);
+ }
+ void run()
+ {
+ tbb::parallel_for(tbb::blocked_range<IndexInt>(0, size), *this);
+ }
+ Grid<T> &me;
+ const Grid<S> &other;
+};
+template<class T, class S> struct gridMult : public KernelBase {
+ gridMult(Grid<T> &me, const Grid<S> &other) : KernelBase(&me, 0), me(me), other(other)
+ {
+ runMessage();
+ run();
+ }
+ inline void op(IndexInt idx, Grid<T> &me, const Grid<S> &other) const
+ {
+ me[idx] *= other[idx];
+ }
+ inline Grid<T> &getArg0()
+ {
+ return me;
+ }
+ typedef Grid<T> type0;
+ inline const Grid<S> &getArg1()
+ {
+ return other;
+ }
+ typedef Grid<S> type1;
+ void runMessage()
+ {
+ debMsg("Executing kernel gridMult ", 3);
+ debMsg("Kernel range"
+ << " x " << maxX << " y " << maxY << " z " << minZ << " - " << maxZ << " ",
+ 4);
+ };
+ void operator()(const tbb::blocked_range<IndexInt> &__r) const
+ {
+ for (IndexInt idx = __r.begin(); idx != (IndexInt)__r.end(); idx++)
+ op(idx, me, other);
+ }
+ void run()
+ {
+ tbb::parallel_for(tbb::blocked_range<IndexInt>(0, size), *this);
+ }
+ Grid<T> &me;
+ const Grid<S> &other;
+};
+template<class T, class S> struct gridDiv : public KernelBase {
+ gridDiv(Grid<T> &me, const Grid<S> &other) : KernelBase(&me, 0), me(me), other(other)
+ {
+ runMessage();
+ run();
+ }
+ inline void op(IndexInt idx, Grid<T> &me, const Grid<S> &other) const
+ {
+ me[idx] /= other[idx];
+ }
+ inline Grid<T> &getArg0()
+ {
+ return me;
+ }
+ typedef Grid<T> type0;
+ inline const Grid<S> &getArg1()
+ {
+ return other;
+ }
+ typedef Grid<S> type1;
+ void runMessage()
+ {
+ debMsg("Executing kernel gridDiv ", 3);
+ debMsg("Kernel range"
+ << " x " << maxX << " y " << maxY << " z " << minZ << " - " << maxZ << " ",
+ 4);
+ };
+ void operator()(const tbb::blocked_range<IndexInt> &__r) const
+ {
+ for (IndexInt idx = __r.begin(); idx != (IndexInt)__r.end(); idx++)
+ op(idx, me, other);
+ }
+ void run()
+ {
+ tbb::parallel_for(tbb::blocked_range<IndexInt>(0, size), *this);
+ }
+ Grid<T> &me;
+ const Grid<S> &other;
+};
+template<class T, class S> struct gridAddScalar : public KernelBase {
+ gridAddScalar(Grid<T> &me, const S &other) : KernelBase(&me, 0), me(me), other(other)
+ {
+ runMessage();
+ run();
+ }
+ inline void op(IndexInt idx, Grid<T> &me, const S &other) const
+ {
+ me[idx] += other;
+ }
+ inline Grid<T> &getArg0()
+ {
+ return me;
+ }
+ typedef Grid<T> type0;
+ inline const S &getArg1()
+ {
+ return other;
+ }
+ typedef S type1;
+ void runMessage()
+ {
+ debMsg("Executing kernel gridAddScalar ", 3);
+ debMsg("Kernel range"
+ << " x " << maxX << " y " << maxY << " z " << minZ << " - " << maxZ << " ",
+ 4);
+ };
+ void operator()(const tbb::blocked_range<IndexInt> &__r) const
+ {
+ for (IndexInt idx = __r.begin(); idx != (IndexInt)__r.end(); idx++)
+ op(idx, me, other);
+ }
+ void run()
+ {
+ tbb::parallel_for(tbb::blocked_range<IndexInt>(0, size), *this);
+ }
+ Grid<T> &me;
+ const S &other;
+};
+template<class T, class S> struct gridMultScalar : public KernelBase {
+ gridMultScalar(Grid<T> &me, const S &other) : KernelBase(&me, 0), me(me), other(other)
+ {
+ runMessage();
+ run();
+ }
+ inline void op(IndexInt idx, Grid<T> &me, const S &other) const
+ {
+ me[idx] *= other;
+ }
+ inline Grid<T> &getArg0()
+ {
+ return me;
+ }
+ typedef Grid<T> type0;
+ inline const S &getArg1()
+ {
+ return other;
+ }
+ typedef S type1;
+ void runMessage()
+ {
+ debMsg("Executing kernel gridMultScalar ", 3);
+ debMsg("Kernel range"
+ << " x " << maxX << " y " << maxY << " z " << minZ << " - " << maxZ << " ",
+ 4);
+ };
+ void operator()(const tbb::blocked_range<IndexInt> &__r) const
+ {
+ for (IndexInt idx = __r.begin(); idx != (IndexInt)__r.end(); idx++)
+ op(idx, me, other);
+ }
+ void run()
+ {
+ tbb::parallel_for(tbb::blocked_range<IndexInt>(0, size), *this);
+ }
+ Grid<T> &me;
+ const S &other;
+};
+template<class T, class S> struct gridScaledAdd : public KernelBase {
+ gridScaledAdd(Grid<T> &me, const Grid<T> &other, const S &factor)
+ : KernelBase(&me, 0), me(me), other(other), factor(factor)
+ {
+ runMessage();
+ run();
+ }
+ inline void op(IndexInt idx, Grid<T> &me, const Grid<T> &other, const S &factor) const
+ {
+ me[idx] += factor * other[idx];
+ }
+ inline Grid<T> &getArg0()
+ {
+ return me;
+ }
+ typedef Grid<T> type0;
+ inline const Grid<T> &getArg1()
+ {
+ return other;
+ }
+ typedef Grid<T> type1;
+ inline const S &getArg2()
+ {
+ return factor;
+ }
+ typedef S type2;
+ void runMessage()
+ {
+ debMsg("Executing kernel gridScaledAdd ", 3);
+ debMsg("Kernel range"
+ << " x " << maxX << " y " << maxY << " z " << minZ << " - " << maxZ << " ",
+ 4);
+ };
+ void operator()(const tbb::blocked_range<IndexInt> &__r) const
+ {
+ for (IndexInt idx = __r.begin(); idx != (IndexInt)__r.end(); idx++)
+ op(idx, me, other, factor);
+ }
+ void run()
+ {
+ tbb::parallel_for(tbb::blocked_range<IndexInt>(0, size), *this);
+ }
+ Grid<T> &me;
+ const Grid<T> &other;
+ const S &factor;
+};
+
+template<class T> struct gridSetConst : public KernelBase {
+ gridSetConst(Grid<T> &grid, T value) : KernelBase(&grid, 0), grid(grid), value(value)
+ {
+ runMessage();
+ run();
+ }
+ inline void op(IndexInt idx, Grid<T> &grid, T value) const
+ {
+ grid[idx] = value;
+ }
+ inline Grid<T> &getArg0()
+ {
+ return grid;
+ }
+ typedef Grid<T> type0;
+ inline T &getArg1()
+ {
+ return value;
+ }
+ typedef T type1;
+ void runMessage()
+ {
+ debMsg("Executing kernel gridSetConst ", 3);
+ debMsg("Kernel range"
+ << " x " << maxX << " y " << maxY << " z " << minZ << " - " << maxZ << " ",
+ 4);
+ };
+ void operator()(const tbb::blocked_range<IndexInt> &__r) const
+ {
+ for (IndexInt idx = __r.begin(); idx != (IndexInt)__r.end(); idx++)
+ op(idx, grid, value);
+ }
+ void run()
+ {
+ tbb::parallel_for(tbb::blocked_range<IndexInt>(0, size), *this);
+ }
+ Grid<T> &grid;
+ T value;
+};
+
+template<class T> template<class S> Grid<T> &Grid<T>::operator+=(const Grid<S> &a)
+{
+ gridAdd<T, S>(*this, a);
+ return *this;
+}
+template<class T> template<class S> Grid<T> &Grid<T>::operator+=(const S &a)
+{
+ gridAddScalar<T, S>(*this, a);
+ return *this;
+}
+template<class T> template<class S> Grid<T> &Grid<T>::operator-=(const Grid<S> &a)
+{
+ gridSub<T, S>(*this, a);
+ return *this;
+}
+template<class T> template<class S> Grid<T> &Grid<T>::operator-=(const S &a)
+{
+ gridAddScalar<T, S>(*this, -a);
+ return *this;
+}
+template<class T> template<class S> Grid<T> &Grid<T>::operator*=(const Grid<S> &a)
+{
+ gridMult<T, S>(*this, a);
+ return *this;
+}
+template<class T> template<class S> Grid<T> &Grid<T>::operator*=(const S &a)
+{
+ gridMultScalar<T, S>(*this, a);
+ return *this;
+}
+template<class T> template<class S> Grid<T> &Grid<T>::operator/=(const Grid<S> &a)
+{
+ gridDiv<T, S>(*this, a);
+ return *this;
+}
+template<class T> template<class S> Grid<T> &Grid<T>::operator/=(const S &a)
+{
+ S rez((S)1.0 / a);
+ gridMultScalar<T, S>(*this, rez);
+ return *this;
+}
+
+//******************************************************************************
+// Other helper functions
+
+// compute gradient of a scalar grid
+inline Vec3 getGradient(const Grid<Real> &data, int i, int j, int k)
+{
+ Vec3 v;
+
+ if (i > data.getSizeX() - 2)
+ i = data.getSizeX() - 2;
+ if (j > data.getSizeY() - 2)
+ j = data.getSizeY() - 2;
+ if (i < 1)
+ i = 1;
+ if (j < 1)
+ j = 1;
+ v = Vec3(data(i + 1, j, k) - data(i - 1, j, k), data(i, j + 1, k) - data(i, j - 1, k), 0.);
+
+ if (data.is3D()) {
+ if (k > data.getSizeZ() - 2)
+ k = data.getSizeZ() - 2;
+ if (k < 1)
+ k = 1;
+ v[2] = data(i, j, k + 1) - data(i, j, k - 1);
+ }
+
+ return v;
+}
+
+// interpolate grid from one size to another size
+
+template<class S> struct knInterpolateGridTempl : public KernelBase {
+ knInterpolateGridTempl(Grid<S> &target,
+ const Grid<S> &source,
+ const Vec3 &sourceFactor,
+ Vec3 offset,
+ int orderSpace = 1)
+ : KernelBase(&target, 0),
+ target(target),
+ source(source),
+ sourceFactor(sourceFactor),
+ offset(offset),
+ orderSpace(orderSpace)
+ {
+ runMessage();
+ run();
+ }
+ inline void op(int i,
+ int j,
+ int k,
+ Grid<S> &target,
+ const Grid<S> &source,
+ const Vec3 &sourceFactor,
+ Vec3 offset,
+ int orderSpace = 1) const
+ {
+ Vec3 pos = Vec3(i, j, k) * sourceFactor + offset;
+ if (!source.is3D())
+ pos[2] = 0; // allow 2d -> 3d
+ target(i, j, k) = source.getInterpolatedHi(pos, orderSpace);
+ }
+ inline Grid<S> &getArg0()
+ {
+ return target;
+ }
+ typedef Grid<S> type0;
+ inline const Grid<S> &getArg1()
+ {
+ return source;
+ }
+ typedef Grid<S> type1;
+ inline const Vec3 &getArg2()
+ {
+ return sourceFactor;
+ }
+ typedef Vec3 type2;
+ inline Vec3 &getArg3()
+ {
+ return offset;
+ }
+ typedef Vec3 type3;
+ inline int &getArg4()
+ {
+ return orderSpace;
+ }
+ typedef int type4;
+ void runMessage()
+ {
+ debMsg("Executing kernel knInterpolateGridTempl ", 3);
+ debMsg("Kernel range"
+ << " x " << maxX << " y " << maxY << " z " << minZ << " - " << maxZ << " ",
+ 4);
+ };
+ void operator()(const tbb::blocked_range<IndexInt> &__r) const
+ {
+ const int _maxX = maxX;
+ const int _maxY = maxY;
+ if (maxZ > 1) {
+ for (int k = __r.begin(); k != (int)__r.end(); k++)
+ for (int j = 0; j < _maxY; j++)
+ for (int i = 0; i < _maxX; i++)
+ op(i, j, k, target, source, sourceFactor, offset, orderSpace);
+ }
+ else {
+ const int k = 0;
+ for (int j = __r.begin(); j != (int)__r.end(); j++)
+ for (int i = 0; i < _maxX; i++)
+ op(i, j, k, target, source, sourceFactor, offset, orderSpace);
+ }
+ }
+ void run()
+ {
+ if (maxZ > 1)
+ tbb::parallel_for(tbb::blocked_range<IndexInt>(minZ, maxZ), *this);
+ else
+ tbb::parallel_for(tbb::blocked_range<IndexInt>(0, maxY), *this);
+ }
+ Grid<S> &target;
+ const Grid<S> &source;
+ const Vec3 &sourceFactor;
+ Vec3 offset;
+ int orderSpace;
+};
+// template glue code - choose interpolation based on template arguments
+template<class GRID> void interpolGridTempl(GRID &target, GRID &source)
+{
+ errMsg("interpolGridTempl - Only valid for specific instantiations");
+}
+
+} // namespace Manta
+#endif