Skip to main content

Faking Buoyancy

For the intro to a game I (perodically) work on, I decided to put the player on a pirate ship. Scrolling textures on a big plane looked pretty boring. To make things a bit angrier I decided to add waves.


The ocean shader displaces the height at each point on a subdivided plane.

float wave(vec3 p) {
	float time = TIME * wave_speed;
	vec2 uv = (p.xz + time * vec2(0.5, 0.5)) * wave_size;
	float dist = length(uv);
	return pow(2.1231, sin(dist) + sin((uv.y + 1.0))) * height;

void vertex() {
	world_position = (MODEL_MATRIX * vec4(VERTEX, 1.0)).xyz;
	VERTEX.y = wave(world_position);

Sum of Sines is a fancy method to simulate water by summing increasingly small amplitude sine waves with different parameters. This is not that, but the idea that summing two sines makes for an interesting wave still holds. I had to play with it a bit until I found something I liked. I simple but interesting wave because that fits the cartoon style.


Next, the boat must stay above the surface of the water. Simply setting the position.y = wave(position) would look odd. Both the position and orientation are affected here.

Starting with the orientation, we can look at this very similarly to looking for the normal of a heightmap. The waves are just a vertex-displaced plane, where the y-displacement is defined by some function. There are a few options for computing the normals:

  • Using small finite differences to approximate the normals at a single point.
  • Use central differences according to the size of the object to calculate the average normals.
  • Using the triangles that make up a quad representing the objects size.

These methods are all very smiple and avoid:

  • Real buoyancy involves volume and dipslacement. This ends up being very noisy and it’s a lot of work.
  • Proper plane fitting algorithms such as Least Squares or Principal Component Analysis. Again the added fidelity (and noise) isn’t worth the effort.

Finite Differences

Finite differences is a method of estimating the slope of a function at a given point. It’s the “rise over run” concept from elementary school. Using a very small step size results in an accurate-enough slope.

$$\Delta x ~= \frac{w(x, z) - w(x + S, z)}{S}$$ $$\Delta z ~= \frac{w(x, z) - w(x, z + S)}{S}$$

To get a normal from these two 2D slopes, we can take the cross product:

$$ \vec n = \begin{bmatrix} S \ \Delta x \ 0 \end{bmatrix} \times \begin{bmatrix} 0 \ \Delta z \ S \end{bmatrix} $$


While the small step gives a very accurate approximation, it might be too accurate. The boat spans many points on the 3D curve. A lot of the boat will clip through the water on concave sections of the curve.

Object-Sized Finite Differences

To mitigate this, we can take the size of our object into account when computing the differences:

$$\Delta x ~= \frac{w(x - 0.5W, z) - w(x + 0.5W, z)}{W}$$ $$\Delta z ~= \frac{w(x, z - 0.5L) - w(x, z + 0.5L)}{L}$$

$$ \vec n = \begin{bmatrix} W \ \Delta x \ 0 \end{bmatrix} \times \begin{bmatrix} 0 \ \Delta z \ L \end{bmatrix} $$

Notice that there are now four samples, and they all go in different directions. This is a version of finite differences called central differences. From these samples, if we rotate our object then have a plane that could fit on the wave but it sits underneath, tangent to the wave’s curve. This is fine in convex scenarios, but it will lie underneath the curve in a concave area. The adjusted position is the mean of our samples.

$$ y’ = \frac{w(x - 0.5W, z) + w(x + 0.5W, z) + w(x, z + .5L) + w(x, z - .5L)}{4} $$

In the case where our rectangle spans across a convex area of the curve, this will move us further down into the water. Using the higher of our adjusted value and the height at the center of the rectangle easily mitigates this.

$$ y’’ = \max(y’, w(x, y)) $$

All of this, in code, looks like this:

func fit_plane(plane: Node3D, size: Vector2) -> Transform3D:
    # take samples
    var center = plane.global_position
	var left = _wave(center + (Vector3.LEFT * size.x / 2).rotated(Vector3.UP, plane.global_rotation.y))
	var right = _wave(center + (Vector3.RIGHT * size.x / 2).rotated(Vector3.UP, plane.global_rotation.y))
	var front = _wave(center + (Vector3.FORWARD * size.y / 2).rotated(Vector3.UP, plane.global_rotation.y))
	var back = _wave(center + (Vector3.BACK * size.y / 2).rotated(Vector3.UP, plane.global_rotation.y))

    # compute the normal
    var dx = right -left
    var dz = back-front
    var normal = -Vector3(size.x, dx, 0).cross(Vector3(0, dz, size.y)).normalized()

    # place the object
    var surface_point = center
    var sample_mean = (left + right + front + back) / 4)
    var surface_point.y = max(_wave(center), sample_mean) 

    # create a transform
	var rotation_axis = Vector3.UP.cross(normal).normalized()
	var rotation_angle = Vector3.UP.angle_to(normal)
	if rotation_axis.length_squared() < .1:
		rotation_axis = Vector3.RIGHT
	return Transform3D(Basis(rotation_axis.normalized(), rotation_angle), surface_point)

Triangles in a Quad

While the central differences approach works okay, it doesn’t handle the subtle rotation across the shorter part of the rectangle. The samples all lie on the midpoints of the sides of the triangles, creating a diamond shape. Sampling at the corners will more accurately represent a rectangular shape, projecting the entire shape onto the curve. The decision to use this approach could be based whether you prefer clipping at the corners of the shape or along the edges.

comparison of trigngle based or central differences

The normal of a triangle in 3D space is the cross-products of two of its sides. There are four triangles formed by the center and any two corners of the quad. The average of just two of the normals yields good-enough results.

In code, this looks like:

func fit_plane(plane: Node3D, Vector3, size: Vector2) -> Transform3D:
	size *= Math.vec2(plane.scale)
    # corners of our rectangle projected onto the plane
    var center = plane.global_position
    var front_r = _wave(center + Vector3(size.x, 0, size.y).rotated(Vector3.UP, plane.global_rotation.y))
    var front_l = _wave(center + Vector3(-size.x, 0, size.y).rotated(Vector3.UP, plane.global_rotation.y))
    var back_r = _wave(center + Vector3(size.x, 0, -size.y).rotated(Vector3.UP, plane.global_rotation.y))
    var back_l = _wave(center + Vector3(-size.x, 0, -size.y).rotated(Vector3.UP, plane.global_rotation.y))

	# front normal
	var v1 = front_l - center
	var v2 = front_r - center
	var normal_f = v1.cross(v2).normalized()

	# back normal
	v1 = back_r - center
	v2 = back_l - center
	var normal_b = v1.cross(v2).normalized()

	# rotation based on average of cross products
	var normal = ((normal_b + normal_f) / 2.0).rotated(Vector3.UP, -plane.global_rotation.y).normalized()

    # place the object
    var surface_point = center
    var sample_mean = (front_r + front_b + back_r + back_l) / 4)
    var surface_point.y = max(_wave(center), sample_mean) 

    # create a transform

The left side uses only normal_b and the right side uses normal_f. The center is the averaged normal.

Riding the Waves

Just for fun, I took things a step further and actually used the waves to drive the motion of the boat. While this may or may for the gameplay to move the boat around, it could be a cool effect for objects that fall into the water.

While I could find the partial derivatives of our wave analytically, I’d have to deal with that each time I change the structure of it. The normals and the gradient are closely related. In fact, when using finite differences, we’ve already computed the gradient.

func fit_plane(plane: Node3D, size: Vector2) -> Transform3D:
    # ...
	var step = .2
	var small_dx = _wave(center + Vector3.RIGHT * step / 2) - _wave(center + Vector3.LEFT * step / 2)
	var small_dz = _wave(center + Vector3.BACK * step / 2) - _wave(center + Vector3.FORWARD * step / 2)
	var grad = (Vector3(small_dx, 0.0, small_dz) / step).normalized()

    var surface_point = center - grad
    # ...

Simply setting the object’s position to the surface_point like this it will slide towards a local minima of the wave and get stuck there.

It also looks a bit too quick. Multiplying by delta time fixes this:

global_position += (surface_transform.origin - global_position) * Vector3(push_strength*delta, 1, push_strength*delta) 

We can add some artistic control by adding a push_strength parameter.

func fit_plane(plane: Node3D, size: Vector2) -> Transform3D:
    # ...
    var surface_point = center - _wave_gradient(center) * strength
    # ...

Because of the smaller step, we don’t keep up with the wave. Because of the wave changing slope at a given point over time, the object can eventually end up changing direction. Instead of getting stuck at some local minimum after encountering one wave, the object looks like it’s swaying back and forth.

Central Differences (Again)

While this is very precise, doing the calculus to get that gradient takes an extra step of manual work. Each time the structure of the _wave function changes, that new function must be differentiated. We can simply re-use the finite differences method from earlier here. The analytic approach could still be useful for validation.

The only modification is how we construct the vector to get the gradient (aka slope) instead of the normal:

var grad = Vector3(dx, 0.0, dz) / step


What about objects that aren’t always in the water? We can re-use all of what we’ve done so far to implement a swimming mechanic. Instead of the model being a child of the plane, we can have a plane that is the child of a body.

First we detect whether or not we’re in the water to turn the influence on or off:

# get swimming position of the parent body
var surface_pos = water.fit_plane(self, Math.vec2(size))
surface_pos.origin -= (global_position - parent.global_position)

# activate swimming mode if we're submerged
if surface_pos.origin.y > global_position.y:
    active = true

# if the player jumps or otherwise ends up above the surface, we're no longer swimming
if parent.global_position.y - surface_pos.origin.y > deactivate_margin:
    active = false

And then apply the influence of the water on the body to velocity rather than directly changing the position.

# align to surface
parent.global_position.y = surface_pos.origin.y

# push the player around with the waves
parent.velocity += (surface_pos.origin - parent.global_position) * delta

This bobbing effect looks neat. After some fine tuning, it could be pretty good, but there is already a lot of motion applied to the player outside their control. Another layer of unpredicatability would take away from the fun, so instead, lets make them stick to the water’s surface.

parent.velocity += (surface_pos.origin - parent.global_position) * Vector3(1, 0, 1) * delta
parent.global_position.y = surface_pos.origin.y

Well, setting the position directly means move_and_slide doesn’t get the opportunity to slide the player and leads to them clipping through walls. Instead, we can assign the y component of velocity so that it puts us exactly on the surface.

var impulse = (surface_pos.origin - parent.global_position) * Vector3(delta, 1/delta, delta) 
parent.velocity += impulse
parent.velocity.y = impulse.y 


I was able to get pretty decent effects without imlementing any hardcore simulations or algorithms. Both the game’s frame-budget and my after work playtime are both very time constrained. Tasks like these where my personal time contraints, gameplay considerations and performance requirements converge are incredibly satisfying.