Palabos程序代码解读(Code Explanation) | tutorial_2 | tutoria_2_2手动设置multiBlock

//部分代码 from tutorial_2_2
template<typename T>
class BounceBackNodes : public DomainFunctional2D {
public:
    BounceBackNodes(plint N, plint radius)
        : cx(N/2),
          cy(N/2),
          innerR(radius),
          outerR(N/2)
    { }
    /// Return true for all cells outside the channel, on which bounce-back
    ///  dynamics must be instantiated.
    virtual bool operator() (plint iX, plint iY) const {
        T rSqr = util::sqr(iX-cx) + util::sqr(iY-cy);
        return rSqr <= innerR*innerR || rSqr >= outerR*outerR;
    }
    virtual BounceBackNodes<T>* clone() const {
        return new BounceBackNodes<T>(*this);
    }
private:
    plint cx;      //< X-position of the center of the half-circle.
    plint cy;      //< Y-position of the center of the half-circle.
    plint innerR;  //< Outer radius of the half-circle.
    plint outerR;  //< Inner radius of the half-circle.
};

图片来源：Palabos Tutorial

这段代码通过virtual bool operator定义了半圆管道以外所有的点，正如下3行所示。

virtual bool operator() (plint iX, plint iY) const {
T rSqr = util::sqr(iX-cx) + util::sqr(iY-cy);
return rSqr <= innerRinnerR || rSqr >= outerRouterR; //innerR= r = N/3

从而将这个域内的cell全部设置为反弹格子。给我们的启发就是可以通过一些线性规划的方程式来设置边界。

不过这个教程文件的重点是下面，即将整个流域划分为4个区域，其中一个区域没有意义以节省计算成本：

图片来源：Palabos Tutorial

//部分代码
// d is the width of the block which is exempted from the full domain.
//d是中间黑色那块域的宽度
plint d = (plint) (2.*std::sqrt((T)util::sqr(radius)-(T)util::sqr(N/4.)));//d=2*(r^2-(N/4)^2)^(1/2)，勾股定理
plint x0 = (N-d)/2 + 1;  // Begin of the exempted block.
plint x1 = (N+d)/2 - 1;  // End of the exempted block.

// Create a block distribution with the three added blocks.
plint envelopeWidth = 1; //“信封”宽度为1，表示网格里的dynamics用的是nearest-neighbor communication类型，意义是这块区域和周围lattice信息交换的格子数
SparseBlockStructure2D sparseBlock(N+1, N/2+1);//整体域，下面通过addBlock来添加里面的域
sparseBlock.addBlock(Box2D(0, x0,      0, N/2),   sparseBlock.nextIncrementalId());//最左的域，id为0，类型为integer
sparseBlock.addBlock(Box2D(x0+1, x1-1, 0, N/4+1), sparseBlock.nextIncrementalId());//中间底下的域，id为1
sparseBlock.addBlock(Box2D(x1, N,      0, N/2),   sparseBlock.nextIncrementalId());//最右的域，id为2

// Instantiate the multi-block, based on the created block distribution and
// on default parameters.
MultiBlockLattice2D<T, DESCRIPTOR> lattice (
    MultiBlockManagement2D (
        sparseBlock,
        defaultMultiBlockPolicy2D().getThreadAttribution(), envelopeWidth ),
    defaultMultiBlockPolicy2D().getBlockCommunicator(),
    defaultMultiBlockPolicy2D().getCombinedStatistics(),
    defaultMultiBlockPolicy2D().getMultiCellAccess<T,DESCRIPTOR>(),
    new BGKdynamics<T,DESCRIPTOR>(omega)
);