Inside boundaries

Question

我正在使用这篇文章来编写流体模拟应用程序。我无法实现内部边界。据我所知，当我为边界内的每个单元格设置边界(在 set_bnd 函数中(时，我应该计算相邻非边界单元格的平均值，如下所示：

for (i = 0 ; i < n ; i++)
{
  for (j = 0 ; j < n ; j++)
  {
     if (isBoundary(i,j)
     {
        sum = 0;
        count = 0;
        if (!isBoundary(i+1,j) {
           sum += x[i+1][j];
        }
        if (!isBoundary(i-1,j) {
           sum += x[i-1][j];
        }
        if (!isBoundary(i,j+1) {
           sum += x[i][j+1];
        }
        if (!isBoundary(i,j-1) {
           sum += x[i-1][j];
        }
        x[i][j] = sum / 4;
     }
  }
}

不幸的是，烟雾被吸收并在与边界表面接触时消失。
我的数学背景不足以理解计算的每个部分，所以如果有人指出我正确的方向，我将不胜感激。

Answer 1

这里有一些代码可以进一步解释。
insideBound是数组(1 - 边界，0 - 空，流体可以通过波谷(

#define FOR_EACH_CELL for ( i=1 ; i<=N ; i++ ) { for ( j=1 ; j<=N ; j++ ) {

void set_bnd ( int N, int b, float * x, int * insideBound )
{
int i, j;
float sum;
int count;
for ( i=1 ; i<=N ; i++ ) {
    x[IX(0  ,i)] = b==1 ? -x[IX(1,i)] : x[IX(1,i)];
    x[IX(N+1,i)] = b==1 ? -x[IX(N,i)] : x[IX(N,i)];
    x[IX(i,0  )] = b==2 ? -x[IX(i,1)] : x[IX(i,1)];
    x[IX(i,N+1)] = b==2 ? -x[IX(i,N)] : x[IX(i,N)];
}
x[IX(0  ,0  )] = 0.5f*(x[IX(1,0  )]+x[IX(0  ,1)]);
x[IX(0  ,N+1)] = 0.5f*(x[IX(1,N+1)]+x[IX(0  ,N)]);
x[IX(N+1,0  )] = 0.5f*(x[IX(N,0  )]+x[IX(N+1,1)]);
x[IX(N+1,N+1)] = 0.5f*(x[IX(N,N+1)]+x[IX(N+1,N)]);
if (!b) return;
FOR_EACH_CELL
    sum = 0.0f;
    count = 0;
    if (insideBound[IX(i,j)] == 1)
    {
        if (insideBound[IX(i-1,j)] != 1)
        {
            count++;
            sum = sum + x[IX(i-1,j)];
        }
        if (insideBound[IX(i+1,j)] != 1)
        {
            count++;
            sum = sum + x[IX(i+1,j)];
        }
        if (insideBound[IX(i,j-1)] != 1)
        {
            count++;
            sum = sum + x[IX(i, j-1)];
        }
        if (insideBound[IX(i,j+1)] != 1)
        {
            count++;
            sum = sum + x[IX(i, j+1)];
        }
        if (count > 0)
        {
            x[IX(i,j)] = -sum / count;
        } else {
            x[IX(i,j)]  = 0;
        }

    }
END_FOR
}

每本书(工作(：
在第一个循环中是设置顶部，右侧，底部和左侧边界单元格。由于对于它们来说，只有一个相邻的单元格未绑定，因此单元格将获得其值。(我不知道为什么它对你来说相反，而 V 的值相同(

在第一个循环之后，设置拐角边界值。在这里，他们从相邻的单元格中获取平均值(我猜因为没有相邻的单元格不是边界的，所以他们使用边界单元格(。

我的，无法正常工作：
if (！b( 返回 - 忽略密度计算并仅更新速度。
循环计算所有边界像元的值(同样，来自本身不是边界的相邻像元的平均值(。我从这种方法中得到了几乎真实的结果，但是密度有很大的损失，还有一些边界太大的错误，流体完全消失了。

Answer 2

我已经设法找到了一个解决方案，这里是为有相同问题的潜在人准备的

void set_bnd ( int N, int b, float * x, int * insideBound )
{
    int i, j;
    float sum, tmp;
    int count;
    for ( i=1 ; i<=N ; i++ ) {
        x[IX(0  ,i)] = b==1 ? -x[IX(1,i)] : x[IX(1,i)];
        x[IX(N+1,i)] = b==1 ? -x[IX(N,i)] : x[IX(N,i)];
        x[IX(i,0  )] = b==2 ? -x[IX(i,1)] : x[IX(i,1)];
        x[IX(i,N+1)] = b==2 ? -x[IX(i,N)] : x[IX(i,N)];
    }
    x[IX(0  ,0  )] = 0.5f*(x[IX(1,0  )]+x[IX(0  ,1)]);
    x[IX(0  ,N+1)] = 0.5f*(x[IX(1,N+1)]+x[IX(0  ,N)]);
    x[IX(N+1,0  )] = 0.5f*(x[IX(N,0  )]+x[IX(N+1,1)]);
    x[IX(N+1,N+1)] = 0.5f*(x[IX(N,N+1)]+x[IX(N+1,N)]);
    if (!b) return;
    for ( i=1 ; i<=N ; i++ ) {
        for ( j=1 ; j<=N ; j++ ) {
            sum = 0.0f;
            count = 0;
            if (insideBound[IX(i,j)] == 1)
            {
                if (insideBound[IX(i-1,j)] != 1)
                {
                    count++;
                    if (b == 2)
                        tmp = -x[IX(i-1,j)];
                    else 
                        tmp = x[IX(i-1,j)];
                    sum = sum + tmp;
                }
                if (insideBound[IX(i+1,j)] != 1)
                {
                    count++;
                    if (b == 2)
                        tmp = -x[IX(i+1,j)];
                    else
                        tmp =  x[IX(i+1,j)];
                    sum = sum + tmp;
                } 
                if (insideBound[IX(i,j-1)] != 1)
                {
                    count++;
                    if (b == 1)
                        tmp = - x[IX(i, j-1)];
                    else 
                        tmp =  x[IX(i, j-1)];
                    sum = sum + tmp;
                } 
                if (insideBound[IX(i,j+1)] != 1)
                {
                    count++;
                    if (b == 1)
                        tmp = -x[IX(i, j+1)];
                    else 
                        tmp = x[IX(i, j+1)];
                    sum = sum + tmp;
                } 
                if (count > 0)
                {
                    x[IX(i,j)] = -sum / count;
                } else {
                    x[IX(i,j)]  = 0;
                }
            }
        }
    }
}

insideBound是布尔数组(0,1(，它表示作为边界的单元格。适用于一个或多个边界区域，但它们应至少为 3 个单元格宽且高度。