1

我正在使用这篇文章来编写一个流体模拟应用程序。我无法实现内部边界。据我所知,当我为边界内的每个单元格设置边界(在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;
     }
  }
}

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

4

2 回答 2

0

我已经设法找到了解决方案,这是为有同样问题的潜在人准备的

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 个单元格。

于 2013-10-15T18:19:55.677 回答
0

这里有一些代码可以进一步解释。
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

}

每本书(工作):
在第一个循环中设置顶部、右侧、底部和左侧边界单元格。因为对于他们来说,只有一个相邻的单元格没有被绑定,所以单元格获得了它的值。(我不知道为什么它对 U 和 V 相同)

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

我的,不能正常工作
如果(!b)返回 - 忽略密度计算并仅更新速度。
循环计算所有边界单元的值(同样,来自不是边界的相邻单元的平均值)。我从这种方法中得到了几乎真实的结果,但是密度损失很大,并且一些边界太大的错误会导致流体完全消失。

于 2013-10-13T11:41:44.357 回答