14,977,910 members
Articles / Multimedia / OpenGL
Article
Posted 23 Apr 2001

319.9K views
135 bookmarked

# Interactive water effect

Rate me:
A short OpenGL demo showing multitexturing effects and clever math optimizations.

## Introduction

This application demonstrates multitexturing effects and clever math optimizations to produce an interactive water effect.

In the calculations used in the code, we work with the following assumptions:

• Incompressible liquid: a given mass of the liquid occupies the same volume all the time whatever the shape it takes. There is conservation of volume.
• The surface of the liquid can be represented by grid heights, i.e. we can see the surface of the liquid as being made up of a series of vertical columns of liquid. We will be able to implement that by an N x M matrix giving the height of each point of surface (N being the number of points in ordinate and M the number of points in x-coordinates). This approximation is valid only for not very deep containers and is tolerable as long as the forces applied to the liquid are not too significant. The liquid cannot splash and a wave cannot break (not bearings).
• The vertical component of the speed of the particles of the liquid can be ignored. The model is not valid any more for stiff slopes or chasms.
• The horizontal component of the speed of the liquid in a vertical column is roughly constant. The movement is uniform within a vertical column. This model is not valid anymore when there are turbulences.

These assumptions enable us to use a simplified form of the equations of the mechanics of fluids. To simplify the reasoning, we first of all will work in 2 dimensions (the height h depends on x-coordinate x and time t):

1. du(x,t)/dt + u(x,t).du(x,t)/dx + g.dh(x,t)/dx = 0
2. 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

where g is gravitational acceleration, h(x, t) is the height of the surface of the liquid, b(x) is the height of the bottom of the container (we suppose that it does not vary according to time t), p(x, t) = h(x, t)-b(x) is the depth of the liquid, and u(x, t) is the horizontal speed of a vertical column of liquid. The equation (1) represents the law of Newton (F=ma) which gives the movement according to the forces applied, and the equation (2) expresses the conservation of volume.

These two equations are nonlinear (nonconstant coefficients) but if the speed of the liquid is small and if the depth varies slowly, we can take the following approximation (we are unaware of the terms multiplied by u(x, t) and p does not depend any more of time t):

3. du(x,t)/dt + g.dh(x,t)/dx = 0
4. dh(x,t)/dt + p(x).du(x,t)/dx = 0

By differentiating the first equation with respect to x and the second equation with respect to t, we obtain:

5. d2u(x,t)/dxdt + g.d2h(x,t)/dx2 = 0
6. d2h(x,t)/dt2 + p(x).d2u(x,t)/dtdx = 0

As u is a "suitable" function (its second derivatives are continuous), its second partial derivatives are equals and we can deduce the following single equation where u is not present:

7. d2h(x,t)/dt2 = gp(x)*d2h(x,t)/dx2

The differential equation of the second member (7) is an equation of wave and expresses the heights' variation as a function of time and the x-coordinate.

In 3 dimensions, the equation (7) has the following form:

8. d2h(x,y,t)/dt2 = gp(x,y) * (d2h(x,y,t)/dx2 + d2h(x,y,t)/dy2)

To be able to work with our heights' grid, we must have a discrete formula, but this one is continuous. Moreover this equation is always nonlinear by the presence of p(x, y). To simplify, we will consider p constant, i.e. a speed of constant wave whatever the depth (we will consider a bottom of flat container, which will limit the "variability" of the depth). We obtain the equation then (9) to discretize:

9. d2h(x,y,t)/dt2 = gp(d2h(x,y,t)/dx2 + d2h(x,y,t)/dy2)

We can discretize the terms of this equation in the following way (using the method of central-differences):

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

By posing Dr = Dx = Dy = resolution of the grid, we obtain the discrete analogue of the equation (9):

10. (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))

where Dt is the variation of time. We can use a more compact notation for this relation of recurrence:

```                                    [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]```

We can then have 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]```

While setting gpDt2/Dr2 = 1/2, we eliminate the last term, and the relation is simplified while giving us:

11. 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) )

This relation of recurrence has a step of 2: the height in t+1 is related to heights in T and in T-1. We can implement that using 2 heights' matrices H1 and 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]

We can then swap the 2 matrices with each step.

To calculate the new height of a point of surface costs only 4 additions, 1 subtraction and a shift-right of 1 (for division by 2).

From the result obtained, we subtract 1/2n of the result (i.e. this same result right shifted of N) to have a certain damping (n=4 gives a rather aesthetic result, n < 4 gives more viscosity, n > 4 gives more fluidity).

Let us notice that the heights of these points are signed, 0 being the neutral level of rest.

From the heights' matrix, we can easily build a polygonal representation by considering each box of the grid as a quadrilateral (or 2 triangles) of which the heights of the 4 vertices are given by the heights of 4 adjacent boxes of the matrix.

In our example, we tessellate our model with `GL_TRIANGLE_STRIP` and we use some multipass effects to get it realistic.

First we perturb a first set of texture coordinates proportionally to vertices normals (the logo's texture coordinates). It looks like refraction.

```/* 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];

}
}
```

Then we calculate a second set of texture coordinates using a fake environment mapping formula (invariant with observer eye's position, just project normals to the xy plane) (the sky's texture coordinates).

```// 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;

}
}```

Then mix the textures together using multitexturing and blending.

```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);```

And to finish the stuff, tessellate using `TRIANGLE_STRIP`s per matrix scanline.

```  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();
}
```

A list of licenses authors might use can be found here

## Share

 Software Developer (Senior) Belgium
No Biography provided

 First PrevNext
 My vote of 5 lovelq52229-Sep-12 5:54 lovelq522 29-Sep-12 5:54
 &%\$#@!\$#%^ nagesh202831-Jul-12 4:22 nagesh2028 31-Jul-12 4:22
 My vote of 5 Mardani Dani5-Jan-12 12:25 Mardani Dani 5-Jan-12 12:25
 need code for iphone dipak.kasabwala6-Apr-11 21:13 dipak.kasabwala 6-Apr-11 21:13
 My vote of 5 BoyDeveloperV230-Oct-10 8:21 BoyDeveloperV2 30-Oct-10 8:21
 Nice BoyDeveloperV230-Oct-10 8:17 BoyDeveloperV2 30-Oct-10 8:17
 Interactive water effect for multitouch. zebulon750183-Jun-10 10:57 zebulon75018 3-Jun-10 10:57
 Re: Interactive water effect for multitouch. Laurent Lardinois3-Jun-10 23:49 Laurent Lardinois 3-Jun-10 23:49
 Great Man .. zebulon7501822-May-10 11:12 zebulon75018 22-May-10 11:12
 iPhone/iPod Touch version on AppStore [modified] Laurent Lardinois27-Mar-10 3:17 Laurent Lardinois 27-Mar-10 3:17
 Re: iPhone/iPod Touch version on AppStore zghotlawala28-Oct-10 17:28 zghotlawala 28-Oct-10 17:28
 did you truly update it for vs2003/vs2005 ? andré_k20-Jun-07 4:06 andré_k 20-Jun-07 4:06
 Re: did you truly update it for vs2003/vs2005 ? Laurent Lardinois8-Feb-10 8:07 Laurent Lardinois 8-Feb-10 8:07
 need your help crazylina1-Feb-07 12:36 crazylina 1-Feb-07 12:36
 Re: A doubt about the use of the code Laurent Lardinois24-Mar-06 8:51 Laurent Lardinois 24-Mar-06 8:51
 Re: A doubt about the use of the code NBee25-Apr-06 11:14 NBee 25-Apr-06 11:14
 cohen Sutherland line clipping algorithm kechalo27-Aug-05 2:31 kechalo 27-Aug-05 2:31
 How to run 2 different(or same) OpenGL objects in one DialogBox? werter13-May-04 21:07 werter1 3-May-04 21:07
 the repartition of light on the ground, under the water walexixis9-Jan-04 17:41 walexixis 9-Jan-04 17:41