交互式水效果。
4.91/5 (53投票s)
2001年4月24日
6分钟阅读
344474
13677
一个简短的OpenGL演示，展示了多纹理效果和巧妙的数学优化。
引言
该应用程序演示了多纹理效果和巧妙的数学优化，以产生交互式水效。
在代码中使用的计算中，我们基于以下假设：
- 不可压缩流体：给定质量的流体在任何形状下都占据相同的体积。体积守恒。
- 流体表面可以用网格高度表示，即我们可以将流体表面看作是由一系列垂直的流体柱组成的。我们可以通过一个NxM矩阵来实现，该矩阵给出表面上每个点的高度（N表示纵坐标上的点数，M表示x坐标上的点数）。这种近似仅适用于不太深的容器，并且只要施加到流体上的力不是非常显著，就可以容忍。流体不会飞溅，波浪也不会破碎（无涌浪）。
- 可以忽略流体粒子的垂直速度分量。该模型对于陡峭的斜坡或深渊不再有效。
- 垂直流体柱中的流体水平速度分量大致恒定。一个垂直流体柱内的运动是均匀的。当存在湍流时，该模型不再有效。
这些假设使我们能够使用简化的流体力学方程。为了简化推理，我们首先在二维（高度h取决于x坐标x和时间t):
- du(x,t)/dt + u(x,t).du(x,t)/dx + g.dh(x,t)/dx = 0
- dp(x,t)/dt + d(u(x,t)p(x,t))/dx = d(h(x,t)-b(x))/dt + p(x,t).du(x,t)/dx + u(x,t).dp(x,t)/dx = 0
其中g是重力加速度，h(x, t)是流体表面的高度，b(x)是容器底的高度（我们假设它不随时间变化t), p(x, t) = h(x, t)-b(x)是流体的深度，而u(x, t)是流体垂直柱的水平速度。方程（1）表示牛顿定律（F=ma），根据作用力给出运动，方程（2）表示体积守恒。
这两个方程都是非线性的（系数非恒定），但如果流体速度很小且深度变化缓慢，我们可以进行以下近似（我们忽略了乘以u(x, t)和p不再依赖于时间 t)
- du(x,t)/dt + g.dh(x,t)/dx = 0
- dh(x,t)/dt + p(x).du(x,t)/dx = 0
通过对第一个方程关于x求导，对第二个方程关于t求导，我们得到
- d2u(x,t)/dxdt + g.d2h(x,t)/dx2 = 0
- d2h(x,t)/dt2 + p(x).d2u(x,t)/dtdx = 0
Asu是一个“合适”的函数（其二阶导数连续），其二阶偏导数相等，我们可以推导出以下单个方程，其中u不存在
- d2h(x,t)/dt2 = gp(x)*d2h(x,t)/dx2
第二个成员（7）的微分方程是一个波动方程，表示高度随时间和x坐标的变化。
在三维中，方程（7）形式如下：
- d2h(x,y,t)/dt2 = gp(x,y) * (d2h(x,y,t)/dx2 + d2h(x,y,t)/dy2)
为了能够使用我们的高度网格，我们必须有一个离散公式，但这是连续的。此外，由于存在p(x, y)，该方程总是非线性的。为简化起见，我们将考虑p常数，即波速不随深度变化（我们将考虑平底容器，这将限制深度的“可变性”）。然后我们得到方程（9）进行离散化
- d2h(x,y,t)/dt2 = gp(d2h(x,y,t)/dx2 + d2h(x,y,t)/dy2)
我们可以使用中心差分法离散化方程的各项：
d2h(x,y,t)/dt2 => (Dh(x,y,t+1) - Dh(x,y,t))/Dt2 = (h(x,y,t+1) - h(x,y,t) - h(x,y,t) + h(x,y,t-1)) / Dt2
d2h(x,y,t)/dx2 => (Dh(x+1,y,t) - Dh(x,y,t))/Dx2 = (h(x+1,y,t) - 2h(x,y,t) + h(x-1,y,t))/Dx2
d2h(x,y,t)/dy2 => (Dh(x,y+1,t) - Dh(x,y,t))/Dy2 = (h(x,y+1,t) - 2h(x,y,t) + h(x,y-1,t))/Dy2
令Dr = Dx = Dy= 网格分辨率，我们得到方程（9）的离散模拟
- (h(x,y,t+1) + h(x,y,t-1) - 2h(x,y,t))/Dt2 = gp/Dr2 . (h(x+1,y,t)+h(x-1,y,t)+h(x,y+1,t)+h(x,y-1,t)-4h(x,y,t))
其中Dt是时间的增量。我们可以为这个递推关系使用更紧凑的记法
[0 1 0] 1/Dt<sup>2</sup> * [1 -2 1] h(x,y) = gp/Dr<sup>2</sup> . [1 -4 1] h(t) [0 1 0]然后我们可以得到h(x, y, t+1):
[0 1 0] h(x,y,t+1) = gpDt<sup>2</sup>/Dr<sup>2</sup> . [1 -4 1]h(t) -h(x,y,t-1) + 2h(x,y,t) [0 1 0] [0 1 0] h(x,y,t+1) = gpDt<sup>2</sup>/Dr<sup>2</sup> . [1 -4 1]h(t) -h(x,y,t-1) + 2h(x,y,t) [0 1 0] [0 1 0] h(x,y,t+1) = gpDt<sup>2</sup>/Dr<sup>2</sup> . [1 0 1]h(t) -h(x,y,t-1) + 2h(x,y,t) - 4gpDt<sup>2</sup>/Dr<sup>2</sup> . h(x,y,t) [0 1 0] [0 1 0] h(x,y,t+1) = gpDt<sup>2</sup>/Dr<sup>2</sup> . [1 0 1]h(t) -h(x,y,t-1) + (2 - 4gpDt<sup>2</sup>/Dr<sup>2</sup>) . h(x,y,t) [0 1 0]通过设置gpDt2/Dr2 = 1/2，我们消除了最后一项，递推关系得以简化，得到
- h(x,y,t+1) = 1/2 * ( h(x+1,y,t) + h(x-1,y,t) + h(x,y+1,t) + h(x,y-1,t) - h(x,y,t-1) )
这个递推关系的步长为2：在t+1的高度与在T和在T-1的高度相关。我们可以使用两个高度矩阵来实现这一点：H1和H2: H2[x, y] = 1/2(H1[x+1, y] + H1[x-1, y] + H1[x, y+1] + H1[x, y-1]) - H2[x, y]
然后我们可以在每一步交换这两个矩阵。
计算一个表面点的新高度只需要4次加法，1次减法和1次右移（除以2）。
从得到的结果中，我们减去1/2n的结果（即该结果右移N）以获得一定的阻尼效果（n=4会得到相当美观的效果，n < 4会产生更大的粘滞性，n > 4会产生更大的流动性）。
请注意，这些点的高度是有符号的，0是静止的平衡水平。
从高度矩阵出发，我们可以很容易地构建一个多边形表示，将网格中的每个盒子视为一个四边形（或2个三角形），其4个顶点的imerick 高度由矩阵中4个相邻盒imerick 的高度给出。
在我们的例子中，我们使用 GL_TRIANGLE_STRIP 对模型进行网格划分，并使用一些多通道效果使其看起来逼真。
首先，我们根据顶点法线（logo的纹理坐标）成比例地扰动一组纹理坐标。这看起来像是折射。
/* calculate ourself normalization */ for(x=0; x < FLOTSIZE*2; x++) { for(y=0; y < FLOTSIZE*2; y++) { sqroot = sqrt(_snormal[x][y][0]*_snormal[x][y][0] + _snormal[x][y][1]*_snormal[x][y][1] + 0.0016f); _snormaln[x][y][0] = _snormal[x][y][0]/sqroot; _snormaln[x][y][1] = _snormal[x][y][1]/sqroot; _snormaln[x][y][2] = 0.04f/sqroot; // perturbate coordinates of background // mapping with the components X,Y of normals... // simulate refraction _newuvmap[x][y][0] = _uvmap[x][y][0] + 0.05f*_snormaln[x][y][0]; _newuvmap[x][y][1] = _uvmap[x][y][1] + 0.05f*_snormaln[x][y][1]; } }
然后，我们使用一个假的反射映射公式（与观察者视线位置无关，只需将法线投影到xy平面）计算第二组纹理坐标（天空的纹理坐标）。
// really simple version of a fake envmap generator for(x=0; x < FLOTSIZE; x++) { for(y=0; y < FLOTSIZE; y++) { // trick : xy projection of normals -> // assume reflection in direction of the normals // looks ok for non-plane surfaces _envmap[x][y][0] = 0.5f + _snormaln[x][y][0]*0.45f; _envmap[x][y][1] = 0.5f + _snormaln[x][y][1]*0.45f; } }
然后使用多纹理和混合来混合纹理。
glEnable(GL_BLEND); // use texture alpha-channel for blending glBlendFunc(GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA); glActiveTextureARB(GL_TEXTURE0_ARB); glBindTexture(GL_TEXTURE_2D, 2); // 2nd texture -> background .. glActiveTextureARB(GL_TEXTURE1_ARB); glBindTexture(GL_TEXTURE_2D, 1); // 2nd texture -> envmap glTexEnvi(GL_TEXTURE_ENV, GL_TEXTURE_ENV_MODE, GL_DECAL); /* enable texture mapping and specify ourself texcoords */ glDisable(GL_TEXTURE_GEN_S); glDisable(GL_TEXTURE_GEN_T);
最后，使用每行扫描线的 TRIANGLE_STRIP 进行网格划分。
for(x=0; x < strip_width; x++) { glBegin(GL_TRIANGLE_STRIP); // WARNING: glTexCoord2fv BEFORE glVertex !!! glMultiTexCoord2fvARB(GL_TEXTURE0_ARB, _newuvmap[x+1][1]); glMultiTexCoord2fvARB(GL_TEXTURE1_ARB, _envmap[x+1][1]); glVertex3fv(_sommet[x+1][1]); // otherwise everything is scrolled !!! for(y=1; y < strip_width; y++) { glMultiTexCoord2fvARB(GL_TEXTURE0_ARB, _newuvmap[x][y]); glMultiTexCoord2fvARB(GL_TEXTURE1_ARB, _envmap[x][y]); glVertex3fv(_sommet[x][y]); glMultiTexCoord2fvARB(GL_TEXTURE0_ARB, _newuvmap[x+1][y+1]); glMultiTexCoord2fvARB(GL_TEXTURE1_ARB, _envmap[x+1][y+1]); glVertex3fv(_sommet[x+1][y+1]); } glMultiTexCoord2fvARB(GL_TEXTURE0_ARB, _newuvmap[x][y]); glMultiTexCoord2fvARB(GL_TEXTURE1_ARB, _envmap[x][y]); glVertex3fv(_sommet[x][y]); glEnd(); }