我浅水实现中出现了奇怪的振荡波纹

7

我一直在尝试在Unity中实现浅水方程,但是遇到了一个奇怪的bug。我的水面出现了奇怪的振荡波纹。我拍了一些截图:

enter image description here enter image description here

这里有一个视频:https://www.youtube.com/watch?v=crXLrvETdjA

我的代码基于Xing Mei的论文Fast Hydraulic Erosion Simulation and Visualization on GPU。你可以在这里找到整个求解器的代码:http://pastebin.com/JktpizHW(或者见下文)。每次我使用论文中的公式时,我都会在注释中添加其编号。

我尝试了不同的时间步长,对于视频,我使用了0.02,将其降低只会使它振荡得更慢。我还尝试了更大的网格(视频使用100,我尝试了200,但水波纹只是更小了)。我多次检查了所有公式,但没有发现任何错误。

有人能找出问题出在哪里吗?

额外信息:

如您在pastebin中所见,我使用C#编写了此程序。我使用Unity作为可视化引擎,并且只是使用网格网格来可视化水。我更改网格的顶点y分量以匹配我计算的高度。
DoUpdate方法获取一个float [][] lowerLayersHeight参数,基本上是水下地形的高度。在视频中,它只是所有0。
public override void DoUpdate(float dt, float dx, float[][] lowerLayersHeight) {
        int x, y;
        float totalHeight, dhL, dhR, dhT, dhB;
        float dt_A_g_l = dt * _A * g / dx; //all constants for equation 2
        float K; // scaling factor for the outflow flux
        float dV;

        for (x=1 ; x <= N ; x++ ) {
                for (y=1 ; y <= N ; y++ ) {
                        //
                        // 3.2.1 Outflow Flux Computation
                        // --------------------------------------------------------------
                        totalHeight = lowerLayersHeight[x][y] + _height[x][y];
                        dhL = totalHeight - lowerLayersHeight[x-1][y] - _height[x-1][y]; //(3)
                        dhR = totalHeight - lowerLayersHeight[x+1][y] - _height[x+1][y];
                        dhT = totalHeight - lowerLayersHeight[x][y+1] - _height[x][y+1];
                        dhB = totalHeight - lowerLayersHeight[x][y-1] - _height[x][y-1];

                        _tempFlux[x][y].left =   Mathf.Max(0.0f, _flux[x][y].left        + dt_A_g_l * dhL ); //(2)
                        _tempFlux[x][y].right =  Mathf.Max(0.0f, _flux[x][y].right       + dt_A_g_l * dhR );
                        _tempFlux[x][y].top =    Mathf.Max(0.0f, _flux[x][y].top         + dt_A_g_l * dhT );
                        _tempFlux[x][y].bottom = Mathf.Max(0.0f, _flux[x][y].bottom      + dt_A_g_l * dhB );

                        float totalFlux = _tempFlux[x][y].left + _tempFlux[x][y].right + _tempFlux[x][y].top + _tempFlux[x][y].bottom;
                        if (totalFlux > 0) {
                                K = Mathf.Min(1.0f, _height[x][y] * dx * dx / totalFlux / dt);  //(4)

                                _tempFlux[x][y].left =   K * _tempFlux[x][y].left;  //(5)
                                _tempFlux[x][y].right =  K * _tempFlux[x][y].right;
                                _tempFlux[x][y].top =    K * _tempFlux[x][y].top;
                                _tempFlux[x][y].bottom = K * _tempFlux[x][y].bottom;
                        }
                        //swap temp and the real one after the for-loops

                        //
                        // 3.2.2 Water Surface
                        // ----------------------------------------------------------------------------------------
                        dV = dt * (
                                //sum in
                                _tempFlux[x-1][y].right + _tempFlux[x][y-1].top + _tempFlux[x+1][y].left + _tempFlux[x][y+1].bottom
                                //minus sum out
                                - _tempFlux[x][y].right - _tempFlux[x][y].top - _tempFlux[x][y].left - _tempFlux[x][y].bottom
                                ); //(6)
                        _tempHeight[x][y] = _height[x][y] + dV / (dx*dx); //(7)
                        //swap temp and the real one after the for-loops
                }
        }
        Helpers.Swap(ref _tempFlux, ref _flux);
        Helpers.Swap(ref _tempHeight, ref _height);
}

Stackoverflow建议您在答案文本中发布最小的代码片段,以重现问题或是最相关的部分。这样可以使问题在未来对谷歌搜索用户保持相关性;如果您链接到pastebin,那么pastebin可能会被删除,这个问题就变得对未来用户毫无价值,因为他们看不到您的代码。 - antiduh
没错,虽然它在算法的某个地方,但我不知道在哪里。它的代码量不是很多,所以我会在这里直接发布它。 - The Oddler
1个回答

6

我自己解决了问题!在开车去朋友家的路上想到了解决办法。问题其实很简单,我在出现错误的代码中对于每一个单元格(或网格点)都计算通量,然后计算高度再转移到下一个单元格。我应该先计算所有单元格的通量,然后在第二次迭代时计算它们的高度。因此代码变为:

for (x=1 ; x <= N ; x++ ) {
    for (y=1 ; y <= N ; y++ ) {
        //
        // 3.2.1 Outflow Flux Computation
        // --------------------------------------------------------------
        ***
    }
}

for (x=1 ; x <= N ; x++ ) {
    for (y=1 ; y <= N ; y++ ) {
        //
        // 3.2.2 Water Surface
        // ---------------------------------------------------------------------------
        ***
    }
}
Helpers.Swap(ref _tempFlux, ref _flux);
Helpers.Swap(ref _tempHeight, ref _height);

当然,***是上述问题中相应的代码。
现在我有一个可用的水模拟。

这是一个非常常见的错误,它应该有自己的名字 :) 例如,在康威的生命游戏中,这是一个经典案例,我在各种数值模拟中也看到过几次。 - Tim Barrass
你好,能否提供更多关于如何运行这个项目的细节?我需要做一个类似的项目,而且我是C#开发人员,我希望能够用C#来完成这个项目,而不是用C++。 - Musaab

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接