Verlet integration 韋爾萊積分 (Part. 1)

Verlet integration:

第一次瞭解到 Verlet integration 像是接觸到了新大陸，原來dynamic sim可以變得如此簡單。

Introduction

When simulating a particle movement, I naturally think of Euler integration, where the position of the particle can be expressed as such: (free falling motion)

$$x=0.5\ast g\ast t^2$$
This is an extremly easy use case. The position of the particle can be represented using a function of time. But in real life, a particle can be affected by multiple forces (and even unstabled force caused non-constant acceleration) and Euler integration also has inaccurate estimation when time step is very large.

In my understanding, the use of Verlet is good for complex particle movement. What it essentially does is relate force with position rather than velocity. The process first calculate the current particle position and the position one step back:

$$velocity=Xcurrent-Xprevious;$$
$$Xprevious=Xcurrent;Xcurrent=Xcurrent+velocity$$

This works magically, here’s some sample code to help understand it.

// without other influence, the following updates the particle position
velocity = posNow - posOld;
posOld 	 = posNow;
posNow  += velocity;

Calculate Force

What if we want to add gravity to the particle like the free falling above? Just add the acceleration to the current position, and the change will add-up on the next integration. The following also applies to forces like friction.

posNow  += GRAVITY; // not sure if we should multipy Time.deltatime

Constraint

Constraint are usually in the form of constraining particles to a specific distance, like spring. When mulitple constraints acts on single particle, one constraint will affect the other. To solve this, loop the constraint multiple time for each time step.

for (int count = 0; count < ITERATION; count++){
		Constriants();
}

// constraint on rope to make fixed distance between segments
void Constriants()
{
	for(int index = 0; index < segment-1; index++)
	{
		float distance = (currentSeg.posNow - nextSeg.posNow).magnitude;
		float error = Mathf.Abs(distance - ropeDist);
		Vector2 changeDir = Vector2.zero;

		if (distance > ropeDist)
			changeDir = (currentSeg.posNow - nextSeg.posNow).normalized;
		else if (distance < ropeDist)
			changeDir = (nextSeg.posNow - currentSeg.posNow).normalized;

		Vector2 changeAmount = changeDir * error;
		if (index == 0)
			nextSeg.posNow += changeAmount;
		else{
			currentSeg.posNow -= changeAmount * 0.5f;
			nextSeg.posNow += changeAmount * 0.5f;
		}
	}
}

Collision

Penalty based system or projection collision reaction

// example for bounding box collision
Vector2 velocity = posNow - posOld;

if (posNow.x > screenWidth){
	posNow.x = screenWidth;
	posOld.x = posNow.x + velocity.x;
}
else if (posNow.x < 0){
	posNow.x = 0;
	posOld.x = posNow.x + velocity.x;
}
if (posNow.y > screenWidth){
	posNow.y = screenWidth;
	posOld.y = posNow.y + velocity.y;
}
else if (posNow.y < 0){
	posNow.y = 0;
	posOld.y = posNow.y + velocity.y;
}

Demo

https://sheldonlei.itch.io/rope-simulation

Resources

https://youtu.be/3HjO_RGIjCU
https://www.youtube.com/watch?v=FcnvwtyxLds
https://en.wikipedia.org/wiki/Verlet_integration#Velocity_Verlet