# Automatic Optimization for Large-Scale Real-Time Coastal Water Simulation

Authors: Shunli Wang; Fengju Kang
Journal: Mathematical Problems in Engineering
 (2016)

---

## Abstract

We introduce an automatic optimization approach for the simulation of large-scale coastal water. To solve the singular problem of water waves obtained with the traditional model, a hybrid deep-shallow-water model is estimated by using an automatic coupling algorithm. It can handle arbitrary water depth and different underwater terrain. As a certain feature of coastal terrain, coastline is detected with the collision detection technology. Then, unnecessary water grid cells are simplified by the automatic simplification algorithm according to the depth. Finally, the model is calculated on Central Processing Unit (CPU) and the simulation is implemented on Graphics Processing Unit (GPU). We show the effectiveness of our method with various results which achieve real-time rendering on consumer-level computer.


---

## Body

## 1. Introduction

Water is a common and important element in nature. It is also frequently encountered in games and other 3D virtual environments. A lot of attention has been paid to rendering water realistically, but the effectiveness of large-scale coastal water is typically limited in games. Large-scale water is often simply modeled to decrease computer cost, so there are not enough details, such as surf, foam, spray, and so forth. More details are needed to ensure the visual effect in today’s computer games and 3D virtual environments.Methods of coastal water simulation can be classified into physically based methods and geometry based methods. Physically based methods mainly solve the shallow-water equations (SWE) with a grid-based or particle-based solver, which requires more than thousands of grid cells or particles [1]. A typical grid-based simulation is developed by Playne et al. [2]. They solve the SWE on multiple graphics processing units (GPUs) and give a standard parallel scheme for realistic fluid simulations. A hybrid grid which consists of regular cells near the water surface and tall cells on the bottom is given by Chentanez and Müller [3]. This method can handle most water phenomena while reducing computation time. It is implemented on a GPU and can get a good real-time performance. Ojeda and Susin [4] simulated shallow-water with the algorithm based on the Lattice Boltzmann Method (LBM). They enhanced the shallow-water scene [5] by adding lower details, surface foam, caustics, and so forth. A fast shallow-water is developed by Abbasov [6]. He constructed a discrete model after numerical simulating a three-dimensional (3D) shallow-water model. This method can greatly simplify the computation by adjusting the cell fill factor upon a sloping shore. The particle-based methods are more suitable for simulating small-scale water with higher details, such as breaking waves [7], water flows into empty regions [8], and interaction of water and rigid objects [9]. Geometry methods usually use a series of shapes to simulate shallow-water. They have been used for a long time. Gerstner model is applied by Fournier and Reeves [10] to the computer graphics field and improved with the elliptic function. However, the water shape he obtained cannot be curved in his simulation. A procedural model for breaking wave in coastal water is provided by De Lima et al. [11]. They generate a wave map and add it to water surface. This method can efficiently simulate the coastal water. Luo et al. [12] put forward a kind of 3D Gerstner wave model with the 3D Bezier curved. These methods can generate a realistic coastal water scene, but they still cannot meet the real-time requirement when large bodies of coastal water need to be rendered. A hybrid shallow-water that combines grid-based and particles based simulation is provided by Chentanez and Müller [13]. Longo et al. [14] gave a review for turbulence in the swash zone and discussed the measurement techniques and generation of breaking waves. There also have shallow-water simulations using a depth-integrated [15], spectrum [16, 17], 2D Boussinesq equations [18], video-based method [19], and so forth. The geometry based methods are suitable for large-scale real-time water because of simple computation. However, these methods cannot simulate the vertical movement, such as curved waves and breaking waves. Although the above methods can get a nice shallow-water scene, they have some limits to meet the requirements of large-scale real-time shallow-water simulations.In order to realize large-scale real-time shallow-water simulations, we analyze the coastal water modeling, height field map generation, ocean wave mesh generation, and rendering. We find some solutions to enhance the real-time performance of the simulations. Our contributions of this paper are as follows:(i)
A construction of the hybrid model that allows the simulation of arbitrary domain scales by coupling of deep water and shallow-water waves.(ii)
An adaptive optimization of water surface grid that can efficiently reduce unnecessary grid cells under the coastal terrain.(iii)
Implementations of our method and water surface effects, such as surface foam and reflection, with the shader language on GPU.This paper is organized as follows. We present our hybrid model in the next section. The optimization of water surface grid is shown in Section3. In Section 4, we present the implementation of our method. Then, we show our rendering results and discussion in Section 5. Finally, we conclude and explore avenues of future work in Section 6.


## 2. Hybrid Coastal Water

When the terrain depth (Depth,D) is greater than half of the wavelength (L), the water is called deep water; otherwise, the water is called shallow water when D is less than L/20, and the intermediate zone shown in Figure 1 is the transition zone.Figure 1
The intermediate zone from deep water to shallow-water.To achieve a realistic coastal water scene, this paper selects the deep water model based on FFT [20] and the shallow-water model based on geometry [21], and the automatic fusion algorithm is used to couple the two water waves in the transitional zone.


### 2.1. Deep Water

The FFT method based on Phillips spectrum first gets a series of double frequency sine waves whose amplitude and phase are relevant and then gets the height of every grid point on the ocean face by calculating with FFT.Assuming that there areM×N discrete points on the xz plane, the vertical offset of the discrete point x=(nLx/N,mLz/M) at time t is(1)hx,t=∑khk,teik·x,where k=kx,ky, kx=2πn/Lx, kz=2πn/Lz, -N/2≤n≤N/2, -M/2≤m≤M/2, and h(k,t) is the height amplitude Fourier and can be calculated with the Phillips spectrum:(2)Phk=A1k4e-1/kL2k^·ω^2.The height amplitude Fourier at timet=0 is(3)h0k=12ξr+iξiPhk,where ξr and ξi are Gaussian random numbers that meet mathematical expectation 0 and variance 1.The height amplitude Fourier at timet is(4)hk,t=h0keiωkt+h0∗-ke-iωkt,where ωk=gk.Actually, in real motion of an ocean wave, the continuous changes of wind direction lead to extrusion between waves which make it become sharp at the top of waves and flat at the bottom following the rolling wave effect. The height map generated by our method is the result of stacking of multiple sine waves. Top of wave is too flat and not real, so we need to add a roll wave model. The model can be realized by changing the position of the grid in the horizontal direction. Here we use the amplitude of the ocean wave to get the offset of discrete pointx in a horizontal direction:(5)Dx,t=-∑kikkhk,teik·x.But this method needs one more FFT, and we simplify it for accelerated computing:(6)Dx,y=SNMx+1,y-Mx,y,Mx,y+1-Mx,y,where S is the size of the FFT grid, N is the number of FFT grid, and M is the wave-height map. The deep water simulation based on FFT is shown in Figure 2(a).Figure 2
The visual results of deep water (a), shallow-water (b), and the hybrid water (c) after coupling.
(a) (b) (c)


### 2.2. Shallow-Water

Compared to the deep water, shallow-water is more complicated due to the effects of seafloor terrain. Its crest line and amplitude will change as the depth decreases. This greatly increases the difficulty of water modeling and rendering. On account of the real-time requirement, this paper uses an improved Fournier model based on geometry instead of the physically based model.Fournier proposed that the motion trails of water points gradually change from circular to elliptical with the decreasing of seabed depth based on the classical Gerstner model. He established the Fournier wave model and applied it to computer graphics. And then he realized the 3D visualization of coastal water. For the sake of rendering of coastal water, including curved waves, Wang et al. [21] made an improvement of Fournier’s model and could quickly simulate coastal water and curved waves at different wind speeds. The formula is as follows:(7)x=x0+RcosαSxsin⁡ϕcos⁡θ+RsinαSzcos⁡ϕcos⁡θ+Δxcosθ,y=y0+RcosαSxsin⁡ϕsin⁡θ+RsinαSzcos⁡ϕsin⁡θ+Δxsinθ,z=z0-RcosαSzcos⁡ϕ+RsinαSxsin⁡ϕ+Δz.Among them, the plane(x,z) is a vertical surface wave, the point (x0,z0) is the static position of a wave point, R is the wave amplitude, ω is the angular frequency, Sx=1/1-e-Kxh and Sz=Sx1-e-Kzh are the major and minor axes of the ellipse, sin⁡α=sin⁡βe-K0h, β is the slope of seafloor terrain, K∞ is the deep water wave number, K(x)=K∞/tanh⁡(K∞h) denotes the wave number at the point (x0,z0), ϕ=-ωt+∑0x0KxΔx, and θ is the wind direction.In order to realize the curved waves, the wave parameters when moving need to be adjusted:Whenϕmin<ϕ≤ϕmin+π(8)ϕ′=Kϕϕ,Δx=Kϕx0,Δz=0,Sx′=121-ϕ-ϕmin-ππ2Sx.Whenϕmin+π<ϕ≤ϕmin+2π(9)ϕ′=ϕ,Δx=ϕ-ϕmin-1.5π0.5π2xcrest-xcurve,Δz=ϕ-ϕmin-1.5π0.5π2zcrest-zcurve,Sx′=Sx.And the limit on the long axis is(10)Sx≤fx,where fx=1/RKxcos⁡α and η is the thin coefficient of the waveform. Its value is 1 by default. The higher the coefficient is, the thinner the waveform is. ϕmin=-tan-1sinαSx/cosαSz denotes the angular ϕ at the wave trough and Kϕ is the curved coefficient with a minimum of 1.The final coastal water is shown in Figure2(b).


### 2.3. Coupling of Water Waves

In order to make the transition of deep water waves and the shallow ones smooth, we use the automatic coupling algorithm to couple the two waves when the depth is in the range of0.05L to 0.5L. The fusion algorithm is shown as follows:(11)hx,y=hFFTx,y,ifd≥0.5Ld-0.05L0.5L-0.05LhFFTx,y+0.5L-d0.5L-0.05Lhgeox,y,if0.05L<d<0.5Lhgeox,y,ifd≤0.05L,where d is the terrain depth of the grid vertex at the current time. hFFT(x,y) and hgeo(x,y) are the heights generated by the deep water and shallow-water model.From the formula, it can be seen that the coupling algorithm will guarantee a smooth transition from deep water to shallow-water. The hybrid water after coupling is shown in Figure2(c).


## 2.1. Deep Water

The FFT method based on Phillips spectrum first gets a series of double frequency sine waves whose amplitude and phase are relevant and then gets the height of every grid point on the ocean face by calculating with FFT.Assuming that there areM×N discrete points on the xz plane, the vertical offset of the discrete point x=(nLx/N,mLz/M) at time t is(1)hx,t=∑khk,teik·x,where k=kx,ky, kx=2πn/Lx, kz=2πn/Lz, -N/2≤n≤N/2, -M/2≤m≤M/2, and h(k,t) is the height amplitude Fourier and can be calculated with the Phillips spectrum:(2)Phk=A1k4e-1/kL2k^·ω^2.The height amplitude Fourier at timet=0 is(3)h0k=12ξr+iξiPhk,where ξr and ξi are Gaussian random numbers that meet mathematical expectation 0 and variance 1.The height amplitude Fourier at timet is(4)hk,t=h0keiωkt+h0∗-ke-iωkt,where ωk=gk.Actually, in real motion of an ocean wave, the continuous changes of wind direction lead to extrusion between waves which make it become sharp at the top of waves and flat at the bottom following the rolling wave effect. The height map generated by our method is the result of stacking of multiple sine waves. Top of wave is too flat and not real, so we need to add a roll wave model. The model can be realized by changing the position of the grid in the horizontal direction. Here we use the amplitude of the ocean wave to get the offset of discrete pointx in a horizontal direction:(5)Dx,t=-∑kikkhk,teik·x.But this method needs one more FFT, and we simplify it for accelerated computing:(6)Dx,y=SNMx+1,y-Mx,y,Mx,y+1-Mx,y,where S is the size of the FFT grid, N is the number of FFT grid, and M is the wave-height map. The deep water simulation based on FFT is shown in Figure 2(a).Figure 2
The visual results of deep water (a), shallow-water (b), and the hybrid water (c) after coupling.
(a) (b) (c)


## 2.2. Shallow-Water

Compared to the deep water, shallow-water is more complicated due to the effects of seafloor terrain. Its crest line and amplitude will change as the depth decreases. This greatly increases the difficulty of water modeling and rendering. On account of the real-time requirement, this paper uses an improved Fournier model based on geometry instead of the physically based model.Fournier proposed that the motion trails of water points gradually change from circular to elliptical with the decreasing of seabed depth based on the classical Gerstner model. He established the Fournier wave model and applied it to computer graphics. And then he realized the 3D visualization of coastal water. For the sake of rendering of coastal water, including curved waves, Wang et al. [21] made an improvement of Fournier’s model and could quickly simulate coastal water and curved waves at different wind speeds. The formula is as follows:(7)x=x0+RcosαSxsin⁡ϕcos⁡θ+RsinαSzcos⁡ϕcos⁡θ+Δxcosθ,y=y0+RcosαSxsin⁡ϕsin⁡θ+RsinαSzcos⁡ϕsin⁡θ+Δxsinθ,z=z0-RcosαSzcos⁡ϕ+RsinαSxsin⁡ϕ+Δz.Among them, the plane(x,z) is a vertical surface wave, the point (x0,z0) is the static position of a wave point, R is the wave amplitude, ω is the angular frequency, Sx=1/1-e-Kxh and Sz=Sx1-e-Kzh are the major and minor axes of the ellipse, sin⁡α=sin⁡βe-K0h, β is the slope of seafloor terrain, K∞ is the deep water wave number, K(x)=K∞/tanh⁡(K∞h) denotes the wave number at the point (x0,z0), ϕ=-ωt+∑0x0KxΔx, and θ is the wind direction.In order to realize the curved waves, the wave parameters when moving need to be adjusted:Whenϕmin<ϕ≤ϕmin+π(8)ϕ′=Kϕϕ,Δx=Kϕx0,Δz=0,Sx′=121-ϕ-ϕmin-ππ2Sx.Whenϕmin+π<ϕ≤ϕmin+2π(9)ϕ′=ϕ,Δx=ϕ-ϕmin-1.5π0.5π2xcrest-xcurve,Δz=ϕ-ϕmin-1.5π0.5π2zcrest-zcurve,Sx′=Sx.And the limit on the long axis is(10)Sx≤fx,where fx=1/RKxcos⁡α and η is the thin coefficient of the waveform. Its value is 1 by default. The higher the coefficient is, the thinner the waveform is. ϕmin=-tan-1sinαSx/cosαSz denotes the angular ϕ at the wave trough and Kϕ is the curved coefficient with a minimum of 1.The final coastal water is shown in Figure2(b).


## 2.3. Coupling of Water Waves

In order to make the transition of deep water waves and the shallow ones smooth, we use the automatic coupling algorithm to couple the two waves when the depth is in the range of0.05L to 0.5L. The fusion algorithm is shown as follows:(11)hx,y=hFFTx,y,ifd≥0.5Ld-0.05L0.5L-0.05LhFFTx,y+0.5L-d0.5L-0.05Lhgeox,y,if0.05L<d<0.5Lhgeox,y,ifd≤0.05L,where d is the terrain depth of the grid vertex at the current time. hFFT(x,y) and hgeo(x,y) are the heights generated by the deep water and shallow-water model.From the formula, it can be seen that the coupling algorithm will guarantee a smooth transition from deep water to shallow-water. The hybrid water after coupling is shown in Figure2(c).


## 3. Water Surface Grid Automatic Simplification

When the water is rendered, grid model is usually used to simulate water surface. Yin et al. [22] divided the surface into a dynamic one and a static one. The dynamic surface is formed by M∗N grid points and revolves around the viewpoint, and the outer space of the dynamic surface is the static surface. Hinsinger et al. [23] presented an adaptive surface grid model. Cui et al. [24] improved the adaptive grid model and put forward a trapezoidal grid model. Johanson and Lejdfors [25] proposed a projection grid model and used it to simulate the water surface in world space.The above studies of wave multiresolution grid did not take the influence of terrain into account. However, shallow-water will not appear under the land terrain, so the surface grid under the land should not be computed and rendered. In order to prevent the unnecessary expense, we proposed an improved multiresolution grid by detecting whether grid points are under the land. Then, rendering of coastal water will be accelerated after reducing computation consuming.


### 3.1. Multiresolution Grid

Multiresolution grid refers to grids with different levels of details (LOD) existing in different regions. It can greatly improve the efficiency of the grid generation without losing the details. Therefore, multiresolution grids are commonly used in real-time rendering of large-scales of water and terrain.In the coastal scene rendering, a part of the water surface grid is always under the terrain. But the former researches usually use the same grid; they did not make changes to the grid. The general grid is shown in Figure3.Figure 3
The general multiresolution grid (a) and the rendering mesh (b).
(a) (b)The grid points covered by terrain will be calculated, generated, and automatically culled by the rendering system. There is no problem in vision, but the process of computing, generating, and culling can severely affect the real-time performance.


### 3.2. Grid Improvement

In order to avoid unnecessary computation, the grid points are detected in advance to ensure whether they are covered by terrain or not when the water surface grid is generating. If they are, the next computation step will be not processed; if not, then continue. The detection will last until the whole grid is generated. Because the coastal terrain is static, the detection does not need to be repeated. Then, the speed of generation and rendering of water can be improved.In the case that coastline data are known, the cover situation can be judged according to the relationship between the grid points and the coastline. The detection process is as follows:Ifxmesh<xline & ymesh<yline, the grid point is not covered by terrain, or it is covered and will not need to be computed.However, the coastline data cannot be generally obtained. Then, the collision detection technology will be used to detect whether the grid point is covered by terrain. The collision detection process is shown in Figure4.Figure 4
The collision detection process.In the process, the collision detection system sends a vertical line from above the grid, and we can manage to keep an account of the height of the point relative to the horizontal line when the line collision with the terrain.Ifzmesh<H, the grid point is not covered by terrain, or it is covered and will not need to be computed. The final grid is given in Figure 5.Figure 5
The improved grid (a) and the rendering result (b).
(a) (b)


## 3.1. Multiresolution Grid

Multiresolution grid refers to grids with different levels of details (LOD) existing in different regions. It can greatly improve the efficiency of the grid generation without losing the details. Therefore, multiresolution grids are commonly used in real-time rendering of large-scales of water and terrain.In the coastal scene rendering, a part of the water surface grid is always under the terrain. But the former researches usually use the same grid; they did not make changes to the grid. The general grid is shown in Figure3.Figure 3
The general multiresolution grid (a) and the rendering mesh (b).
(a) (b)The grid points covered by terrain will be calculated, generated, and automatically culled by the rendering system. There is no problem in vision, but the process of computing, generating, and culling can severely affect the real-time performance.


## 3.2. Grid Improvement

In order to avoid unnecessary computation, the grid points are detected in advance to ensure whether they are covered by terrain or not when the water surface grid is generating. If they are, the next computation step will be not processed; if not, then continue. The detection will last until the whole grid is generated. Because the coastal terrain is static, the detection does not need to be repeated. Then, the speed of generation and rendering of water can be improved.In the case that coastline data are known, the cover situation can be judged according to the relationship between the grid points and the coastline. The detection process is as follows:Ifxmesh<xline & ymesh<yline, the grid point is not covered by terrain, or it is covered and will not need to be computed.However, the coastline data cannot be generally obtained. Then, the collision detection technology will be used to detect whether the grid point is covered by terrain. The collision detection process is shown in Figure4.Figure 4
The collision detection process.In the process, the collision detection system sends a vertical line from above the grid, and we can manage to keep an account of the height of the point relative to the horizontal line when the line collision with the terrain.Ifzmesh<H, the grid point is not covered by terrain, or it is covered and will not need to be computed. The final grid is given in Figure 5.Figure 5
The improved grid (a) and the rendering result (b).
(a) (b)


## 4. GPU Implementation Method

In order to rapidly generate large bodies of near ocean wave scene, graphics hardware shader language is used in rendering. The water height map data are obtained by the calculating of shallow-water model. Then, it is loaded in the vertex shader and is fused with Perlin noise. Final water height and offset in horizontal direction are computed. Each vertex normal vector and color data are calculated in the pixel shader, so the water scene is realized. The surface reflection and foam effect are also added to increase the scene fidelity.


### 4.1. Rendering Based on Height Field Map

The curved model can be computed to obtain the height and horizontal offset of each vertex. The result will be saved in a height map. Each pixel of the height map can store a vertex information. It is expressed by the formula(12)colordata·rgb=x,y,z.The height field map data can be read and interpolated in the vertex shader, and Perlin noise is added; then, the final vertex information is expressed by(13)x,y,z=colordata·r,colordata·g,colordata·b+Noiseperlin.In the pixel shader, the vertex normal vector and color are calculated, the vertex transparency is computed according to the depth:(14)colorfinal·rgb=fresnel∗colorwater·rgb+1-fresnel∗colordifuse,colorfinal·a=1-0.6e-h.Here fresnel denotes the fresnel reflection coefficient.The water surface will be getting more and more transparent with the decreasing of depth. The value of 0 means fully transparent. In this paper, the minimum transparency value is 0.6.


### 4.2. Realization of the Reflection and Foam Effect

The rendering scene above sea level is rendered to a texture through the RTT (render to texture) technology. The reflection texture is loaded in the pixel shader, and then the reflection effect of the wave surface is realized by using dynamic texture mapping technology. The specific calculation formula is as follows:(15)colorfinal·rgb=fresnel∗colorwater·rgb+1-fresnel∗colorreflect.Foam will be intermittently generated at the wave crest when the wave is moving. And foams appear in the trough and crest of the waves during wave run-up and run-down [14]. So we applied two different methods to achieve the effect of foam. We load a transparent foam texture and simulate the dynamic foam effect by controlling the offset speed of texture with texture mapping technology. The implementation detail and rendering result are shown in Box 1 and Figure 6:(16)colorfinal=density∗colorfinal+1-density∗colorfoam.Box 1:The implementation detail.
/ / Combine waterbody color and reflected color vec4 water_color = mix(body_color, reflected_color, fresnel); vec4 final_color = water_color; float depth = texture2D( DeepMap, vec2(pos_local.x, (pos_local.y))/256.0 ).r * 30.0 - 15.0; // Crest Foam float alphaFactor = texture2D( FoamMap, pos_local.xy/64.0 + g_PerlinMovement / 2.0).r; float foamFactor = alphaFactor * clamp(pos_local.z / 2.0 - 0.1, 0.0, 1.0); final_color = vec4(final_color.xyz, (300.0 - pos_local.x) / 300.0 ); final_color = mix(final_color, vec4(1.0, 1.0, 1.0, 1.0), foamFactor); // Trough Foam alphaFactor = texture2D( FoamMap, vec2(pos_local.x/32.0, pos_local.y/32.0) + g_PerlinMovement / 2.0).r; float foamFactor0 = alphaFactor * clamp( (depth + 3.0) / 5.0, 0.0, 1.0); alphaFactor = texture2D( FoamMap, pos_local.xy/64.0 - g_PerlinMovement / 3.0).r; float foamFactor1 = alphaFactor * clamp( -0.6 - pos_local.z, 0.0, 1.0) * clamp( (depth + 20.0) / 15.0, 0.0, 1.0); final_color = mix(final_color, vec4(1.0, 1.0, 1.0, 1.0), foamFactor0); final_color = mix(final_color, vec4(1.0, 1.0, 1.0, 1.0), foamFactor1); // Fog float fogFactor = computeFogFactor( g_FogDensity, gl_FogFragCoord ); water_color = final_color; // Sun spots float cos_spec = clamp(dot(reflect_vec, g_SunDir), 0.0, 1.0); float sun_spot = pow(cos_spec, g_Shineness); water_color.xyz += g_SunColor * sun_spot; gl_FragColor = water_color;Figure 6
The reflection and foam effects.


## 4.1. Rendering Based on Height Field Map

The curved model can be computed to obtain the height and horizontal offset of each vertex. The result will be saved in a height map. Each pixel of the height map can store a vertex information. It is expressed by the formula(12)colordata·rgb=x,y,z.The height field map data can be read and interpolated in the vertex shader, and Perlin noise is added; then, the final vertex information is expressed by(13)x,y,z=colordata·r,colordata·g,colordata·b+Noiseperlin.In the pixel shader, the vertex normal vector and color are calculated, the vertex transparency is computed according to the depth:(14)colorfinal·rgb=fresnel∗colorwater·rgb+1-fresnel∗colordifuse,colorfinal·a=1-0.6e-h.Here fresnel denotes the fresnel reflection coefficient.The water surface will be getting more and more transparent with the decreasing of depth. The value of 0 means fully transparent. In this paper, the minimum transparency value is 0.6.


## 4.2. Realization of the Reflection and Foam Effect

The rendering scene above sea level is rendered to a texture through the RTT (render to texture) technology. The reflection texture is loaded in the pixel shader, and then the reflection effect of the wave surface is realized by using dynamic texture mapping technology. The specific calculation formula is as follows:(15)colorfinal·rgb=fresnel∗colorwater·rgb+1-fresnel∗colorreflect.Foam will be intermittently generated at the wave crest when the wave is moving. And foams appear in the trough and crest of the waves during wave run-up and run-down [14]. So we applied two different methods to achieve the effect of foam. We load a transparent foam texture and simulate the dynamic foam effect by controlling the offset speed of texture with texture mapping technology. The implementation detail and rendering result are shown in Box 1 and Figure 6:(16)colorfinal=density∗colorfinal+1-density∗colorfoam.Box 1:The implementation detail.
/ / Combine waterbody color and reflected color vec4 water_color = mix(body_color, reflected_color, fresnel); vec4 final_color = water_color; float depth = texture2D( DeepMap, vec2(pos_local.x, (pos_local.y))/256.0 ).r * 30.0 - 15.0; // Crest Foam float alphaFactor = texture2D( FoamMap, pos_local.xy/64.0 + g_PerlinMovement / 2.0).r; float foamFactor = alphaFactor * clamp(pos_local.z / 2.0 - 0.1, 0.0, 1.0); final_color = vec4(final_color.xyz, (300.0 - pos_local.x) / 300.0 ); final_color = mix(final_color, vec4(1.0, 1.0, 1.0, 1.0), foamFactor); // Trough Foam alphaFactor = texture2D( FoamMap, vec2(pos_local.x/32.0, pos_local.y/32.0) + g_PerlinMovement / 2.0).r; float foamFactor0 = alphaFactor * clamp( (depth + 3.0) / 5.0, 0.0, 1.0); alphaFactor = texture2D( FoamMap, pos_local.xy/64.0 - g_PerlinMovement / 3.0).r; float foamFactor1 = alphaFactor * clamp( -0.6 - pos_local.z, 0.0, 1.0) * clamp( (depth + 20.0) / 15.0, 0.0, 1.0); final_color = mix(final_color, vec4(1.0, 1.0, 1.0, 1.0), foamFactor0); final_color = mix(final_color, vec4(1.0, 1.0, 1.0, 1.0), foamFactor1); // Fog float fogFactor = computeFogFactor( g_FogDensity, gl_FogFragCoord ); water_color = final_color; // Sun spots float cos_spec = clamp(dot(reflect_vec, g_SunDir), 0.0, 1.0); float sun_spot = pow(cos_spec, g_Shineness); water_color.xyz += g_SunColor * sun_spot; gl_FragColor = water_color;Figure 6
The reflection and foam effects.


## 5. Results and Discussion

To evaluate our method, we have created several typical coastal water scenes on a 2.66 GHz Core i7 and NVIDIA GTX 560. The hybrid water model is computed on CPU, and the rendering is implemented on GPU. The results are shown in Figure7. In Figures 7(a) and 7(b), the thin coefficient is 1.0, while its value is 1.2 in Figures 7(c) and 7(d). Compared with 1.0, the waveform is thinner with the coefficient of 1.2. But there is a problem that the waveforms will be crossed when the coefficient getting higher. It is the limit of our method, and we will research it in future work. The foam effects are different in the crest and trough of waves. In our method, foam at the crest moves forward and foam at the trough moves backward. So we can simulate the huge amount of bubbles fast. Table 1 lists the performance values of different examples. Since the model is geometric and the rendering is implemented with shader language, our method makes it possible to simulate a large-scale coastal water scene at a stable frame rate of 60 frames per second.Table 1
Timing per frame for various examples in millisecond.
Grid size Computation grid size Model Rendering Total 5122 1282 6.56 5.06 12.06 2562 21.67 13.10 36.85 5122 115.89 18.56 134.58  10242 1282 6.58 5.33 13.83 2562 22.36 14.07 38.02 5122 116.23 21.99 138.69Figure 7
The rendering result of coastal water scene with different thin coefficients.
(a)
η = 1.0 (b)
η = 1.0 (c)
η = 1.2 (d)
η = 1.2The comparison with the references’ methods is shown in Table2. Ojeda’s method [4] is based on physics. It can handle arbitrary underlying terrain and interact with dynamic rigid bodies. They compute the water model both on CPU and on GPU, but only the result on GPU can meet the requirement of real-time performance. Another problem of our method is that it can only get the simulation of small-scale shallow-water. Luo’s method [12] used a Bezier curve, so it is based on geometry. They can get a nice real-time middle-scale coastal water after GPU accelerating. However, they need to design a coastline before rendering. Once the rendering is begun, the coastline cannot be changed. As our method can get a high frame rate, we will not need a GPU based computation of the water model. We can also handle any underlying terrain since the model is computed in real-time. Furthermore, we also consider the deep water when the depth is getting larger. It can also help to make the coastal scene more real.Table 2
A comparison of the difference between Luo’s method, Ojeda’s method, and our method.
Ojeda’s method Luo’s method Our method Computation platform CPU GPU GPU CPU Computation grid size 1282 1282 2502 1282 Water scale Small Small Middle Large Averaged timing per frame 410.88 ms 5.28 ms 16.67 ms 12.06 msOne drawback of our method is that the rendering results lack higher details. The geometric based model is not the same as the physically based model, and it cannot simulate the higher details, such as splash and breaking dam. The other problem is that the waveform will become crossed if we do not limit the thin coefficient. However, our method offers a high rendering speed even on an ordinary graphic hardware. This can allow more consumer-level computers to run the coastal water simulation.


## 6. Conclusions

In this paper, we present an automatic optimization approach for the rendering of large-scale real-time coastal water. A hybrid coastal water model is estimated for large-scale coastal water according to the depth of water. The deep water is also considered in our hybrid coastal water. Unnecessary grid cells are simplified by our automatic simplification algorithm. Finally, we implement our method on GPU by using shader language. The reflection and foam effect are also added to the water surface to enhance the final rendering result.There are several ways to further improve the visual appearance of the coastal water. The extraordinary detailed effects like waterfall, breaking dam, and splashes are needed to be complex modeled by using particles systems. The interaction with rigid bodies is also a challenging work in the future.



---

*Source: 9034649-2016-10-31.xml*