# Soft-blob physics

This is the first of two articles on how I created the glass effect in my latest generative art project Vitreous, in this one we’re exploring how I made the shapes for each of the glass pieces by implementing various generative algorithms and physical simulation techniques.

## Starting points

The first step was to figure out sensible starting points for each blob, these won’t necessarily be the centre of each blob, but where the blob will grow out from. I wanted to choose starting points that would result in evenly-sized blobs, but have an unpredictable distribution. There are a million ways of going about this, but here are the ideas I played with.

### Jitter grid

This method involves making a regular grid and then shifting each of the points in a random pattern.

## Code snippet

```
for (let column = 0; column < columns; column++) {
for (let row = 0; row < rows; row++) {
const angle = Math.random() * Math.PI * 2
let x = jitter * Math.cos(angle) + column * gap
let y = jitter * Math.sin(angle) + row * gap
points.push([x, y])
}
}
```

```
for (let column = 0; column < columns; column++) {
for (let row = 0; row < rows; row++) {
const angle = Math.random() * Math.PI * 2
let x = jitter * Math.cos(angle) + column * gap
let y = jitter * Math.sin(angle) + row * gap
points.push([x, y])
}
}
```

With a high level of jitter you can no longer see the underlying grid, but some points are much closer together than others. This is much easier to see if you turn on ‘closest dot boundary’ (this draws circles around the dots which take up as much space as possible without broaching another dot’s circle). The higher the jitter the lower the packing density, which means our blobs will have to expand a lot to take up the empty space, which will result in very uneven blobs – this isn’t necessarily a bad thing, but it’s not ideal for this use case.

By starting with a hexagonal rather than square grid we can make the dot spread look much more natural with a smaller amount of jitter, and still obscure the underlying grid.

We could take this one step further by employing Lloyd’s algorithm to spread the points out more evenly – this uses Voronoi tesellation (more on that in a bit) to divide the layout into cells, and then moves each point to the centroid of its cell. After several iterations the points reach an equilibrium where they’re fairly evenly spaced out but in an irregular formation.

Lloyd’s algorithm seems like a perfect solution to this problem but I wanted to try one of my favourite algorithms first because it has an extra feature that I thought could come in handy later on…

### Poisson-disk sampling

Poisson-disk sampling aims to pack a series of points tightly while maintaining a minimum distance. I use Kevin Chapelier’s poisson-disk-sampling npm package, which is based on Robert Bridson’s 2007 algorithm.

## Code snippet

```
import { PoissonDiskSampling } from 'poisson-disk-sampling'
const sampling = new PoissonDiskSampling({
shape: [width, height],
minDistance: 5,
maxDistance: 10,
})
const points = sampling.fill()
```

```
import { PoissonDiskSampling } from 'poisson-disk-sampling'
const sampling = new PoissonDiskSampling({
shape: [width, height],
minDistance: 5,
maxDistance: 10,
})
const points = sampling.fill()
```

This algorithm results in a really clean, organic spread of points, but if you turn on the closest dot boundary you’ll see there are some quite large gaps between the clusters of circles.

To have a good starting point for my growth system I want there to be as few gaps as possible, with as consistently-sized boundaries as possible, while still looking natural. It’s a tall order.

I experimented with a few ways to close these gaps up without disrupting the organic distribution too much. I settled on a system where each point would move away from its closest circle-packed neighbour. In areas where they’re packed closely this wouldn’t change the order because they would ping between one another, but it does encourage points to move into the empty spaces.

## Code snippet

```
const getDist = (pointA, pointB) =>
Math.sqrt((pointB.x - pointA.x) ** 2 + (pointB.y - pointA.y) ** 2)
const getClosestDist = (pointA, points) => {
let closestDist = Infinity
points.forEach((pointB) => {
if (pointA === pointB) return
const dist = getDist(pointA, pointB)
if (dist < closestDist) closestDist = dist
})
return closestDist
}
const separatePoints = (points) => {
// get the circle-packed distance for each point
const blobRadii = []
points.forEach((point) => {
const dist = getClosestDist(point, points)
blobRadii.push(dist / 2)
})
points.forEach((pointA, iA) => {
let closestDistance = Infinity
let closestPoint = pointA
points.forEach((pointB, iB) => {
if (pointA === pointB) return
// find the point whose circle-packed edge is closest
const dist = getDist(pointA, pointB) - blobRadii[iA] - blobRadii[iB]
if (dist < closestDistance) {
closestDistance = dist
closestPoint = pointB
}
})
// move the points away from each other
const travel = 5 / getDist(pointA, closestPoint)
pointA.x -= (closestPoint.x - pointA.x) * travel
pointA.y -= (closestPoint.y - pointA.y) * travel
closestPoint.x += (closestPoint.x - pointA.x) * travel
closestPoint.y += (closestPoint.y - pointA.y) * travel
})
}
```

```
const getDist = (pointA, pointB) =>
Math.sqrt((pointB.x - pointA.x) ** 2 + (pointB.y - pointA.y) ** 2)
const getClosestDist = (pointA, points) => {
let closestDist = Infinity
points.forEach((pointB) => {
if (pointA === pointB) return
const dist = getDist(pointA, pointB)
if (dist < closestDist) closestDist = dist
})
return closestDist
}
const separatePoints = (points) => {
// get the circle-packed distance for each point
const blobRadii = []
points.forEach((point) => {
const dist = getClosestDist(point, points)
blobRadii.push(dist / 2)
})
points.forEach((pointA, iA) => {
let closestDistance = Infinity
let closestPoint = pointA
points.forEach((pointB, iB) => {
if (pointA === pointB) return
// find the point whose circle-packed edge is closest
const dist = getDist(pointA, pointB) - blobRadii[iA] - blobRadii[iB]
if (dist < closestDistance) {
closestDistance = dist
closestPoint = pointB
}
})
// move the points away from each other
const travel = 5 / getDist(pointA, closestPoint)
pointA.x -= (closestPoint.x - pointA.x) * travel
pointA.y -= (closestPoint.y - pointA.y) * travel
closestPoint.x += (closestPoint.x - pointA.x) * travel
closestPoint.y += (closestPoint.y - pointA.y) * travel
})
}
```

If this algorithm is left running it will eventually settle into a hexagonal pattern, and as it doesn’t have a bounding box the points will continue to expand outwards, but by running the simulation for a few frames we can find a balance between maintaining disorder and eliminating the unnecessary gaps. Huzzah!

The other advantage of using poisson-disk sampling is that the minimum gap between points can be variable, which allows you to modulate the density of points however you want. I used this feature in Vitreous with octaval noise to subtly vary the sizes of the blobs across the composition while maintaining similar sizes between neighbours.

## Blob formation

Now we’ve perfected the algorithm to choose our starting points, it’s onto using them to make 2D shapes. As with all my technical explorations I start with a technique that’s easy to implement and follow my nose from there.

### Voronoi tesellation

A voronoi diagram partitions a plane into polygons by growing out from starting points (‘seeds’). Each cell contains one seed, and any point on the diagram is within the cell of its nearest seed. That technical definition isn’t very intuitive but they make sense when you see them.

My first attempt at making soft blobs was to use a voronoi diagram and round the edges of each cell using quadratic bezier curves. The easiest way to add gaps between the cells is to scale the shapes from their seed points.

## Code snippet

```
import Voronoi from 'voronoi'
const settings = {
scale: 0.9,
rounding: 0.5,
}
const map = (value, min1, max1, min2, max2) =>
((value - min1) / (max1 - min1)) * (max2 - min2) + min2
const voronoi = new Voronoi()
const diagram = voronoi.compute(points, {
xl: 0,
xr: width,
yt: 0,
yb: height,
})
diagram.cells.forEach((cell) => {
const center = cell.site
c.save()
c.translate(center.x, center.y)
c.scale(settings.scale, settings.scale)
c.beginPath()
// start in the middle of the first edge
const startEdge = cell.halfedges[0]
const startPoint = startEdge.getStartpoint()
const endPoint = startEdge.getEndpoint()
c.moveTo(
map(0.5, 0, 1, startPoint.x, endPoint.x) - center.x,
map(0.5, 0, 1, startPoint.y, endPoint.y) - center.y
)
// loop through edges
cell.halfedges.forEach((edge, edgeI) => {
const nextEdge = cell.halfedges[(edgeI + 1) % cell.halfedges.length]
const startPoint = edge.getStartpoint()
const pivotPoint = edge.getEndpoint()
const endPoint = nextEdge.getEndpoint()
// start quadratic curve between middle of edge and end of edge
const openX = map(settings.rounding, 0, 2, pivotPoint.x, startPoint.x)
const openY = map(settings.rounding, 0, 2, pivotPoint.y, startPoint.y)
// end quadratic curve between end of edge and middle of next edge
const closeX = map(settings.rounding, 0, 2, pivotPoint.x, endPoint.x)
const closeY = map(settings.rounding, 0, 2, pivotPoint.y, endPoint.y)
c.lineTo(openX - center.x, openY - center.y)
c.quadraticCurveTo(
pivotPoint.x - center.x,
pivotPoint.y - center.y,
closeX - center.x,
closeY - center.y
)
})
c.closePath()
c.stroke()
c.restore()
})
```

```
import Voronoi from 'voronoi'
const settings = {
scale: 0.9,
rounding: 0.5,
}
const map = (value, min1, max1, min2, max2) =>
((value - min1) / (max1 - min1)) * (max2 - min2) + min2
const voronoi = new Voronoi()
const diagram = voronoi.compute(points, {
xl: 0,
xr: width,
yt: 0,
yb: height,
})
diagram.cells.forEach((cell) => {
const center = cell.site
c.save()
c.translate(center.x, center.y)
c.scale(settings.scale, settings.scale)
c.beginPath()
// start in the middle of the first edge
const startEdge = cell.halfedges[0]
const startPoint = startEdge.getStartpoint()
const endPoint = startEdge.getEndpoint()
c.moveTo(
map(0.5, 0, 1, startPoint.x, endPoint.x) - center.x,
map(0.5, 0, 1, startPoint.y, endPoint.y) - center.y
)
// loop through edges
cell.halfedges.forEach((edge, edgeI) => {
const nextEdge = cell.halfedges[(edgeI + 1) % cell.halfedges.length]
const startPoint = edge.getStartpoint()
const pivotPoint = edge.getEndpoint()
const endPoint = nextEdge.getEndpoint()
// start quadratic curve between middle of edge and end of edge
const openX = map(settings.rounding, 0, 2, pivotPoint.x, startPoint.x)
const openY = map(settings.rounding, 0, 2, pivotPoint.y, startPoint.y)
// end quadratic curve between end of edge and middle of next edge
const closeX = map(settings.rounding, 0, 2, pivotPoint.x, endPoint.x)
const closeY = map(settings.rounding, 0, 2, pivotPoint.y, endPoint.y)
c.lineTo(openX - center.x, openY - center.y)
c.quadraticCurveTo(
pivotPoint.x - center.x,
pivotPoint.y - center.y,
closeX - center.x,
closeY - center.y
)
})
c.closePath()
c.stroke()
c.restore()
})
```

Although this method produces blobs that look either too rigid or too rounded, I was surprised by how effective it was for a first attempt. The most obvious flaw is on short edges where there isn’t much of a line to curve around so it comes to a point.

My next development was to implement a system that mimics how voronoi diagrams are calculated – treat each seed as an expanding circle that stops when it meets the boundary of another circle. In practical terms I implemented this with a series of particles that move outwards from their seed point and stop when they’re close to a particle from another blob. This allows us to add a more accurate border between blobs by stopping each particle at a certain distance from other particles.

## Code snippet

```
import * as THREE from 'three'
const temp = new THREE.Vector2()
const settings = {
edgePointCount: 20,
stoppingDistance: [5, 40],
}
const map = (value, min1, max1, min2, max2, clampResult) => {
const result = ((value - min1) / (max1 - min1)) * (max2 - min2) + min2
return clampResult ? Math.max(min2, Math.min(result, max2)) : result
}
class Blob {
constructor(x, y) {
this.centre = new THREE.Vector2(x, y)
this.particles = []
for (let i = 0; i < settings.edgePointCount; i++) {
const angle = (i / settings.edgePointCount) * Math.PI * 2
const x = Math.cos(angle) + this.centre.x
const y = Math.sin(angle) + this.centre.y
this.particles.push(new THREE.Vector2(x, y))
}
}
repelBlobs(blobs) {
this.particles.forEach((particle) => {
// find closest particle
let closestDist = Infinity
blobs.forEach((blobB) => {
if (this === blobB) return
// optional: ignore far-away blobs here for performance
blobB.particles.forEach((particleB) => {
const dist = particle.distanceTo(particleB)
if (dist < closestDist) closestDist = dist
})
})
// calculate the vector between particle and centre
temp.set(particle.x, particle.y)
temp.sub(this.centre)
temp.normalise()
temp.multiplyScalar(
map(
closestDist,
settings.stoppingDistance[0],
settings.stoppingDistance[1],
0,
0.2,
true
)
)
particle.add(temp)
})
}
draw() {
blobs.forEach((blob) => {
c.beginPath()
c.moveTo(blob.particles[0].x, blob.particles[0].y)
blob.particles.slice(1).forEach((particle) => {
c.lineTo(particle.x, particle.y)
})
c.closePath()
c.fill()
})
}
}
```

```
import * as THREE from 'three'
const temp = new THREE.Vector2()
const settings = {
edgePointCount: 20,
stoppingDistance: [5, 40],
}
const map = (value, min1, max1, min2, max2, clampResult) => {
const result = ((value - min1) / (max1 - min1)) * (max2 - min2) + min2
return clampResult ? Math.max(min2, Math.min(result, max2)) : result
}
class Blob {
constructor(x, y) {
this.centre = new THREE.Vector2(x, y)
this.particles = []
for (let i = 0; i < settings.edgePointCount; i++) {
const angle = (i / settings.edgePointCount) * Math.PI * 2
const x = Math.cos(angle) + this.centre.x
const y = Math.sin(angle) + this.centre.y
this.particles.push(new THREE.Vector2(x, y))
}
}
repelBlobs(blobs) {
this.particles.forEach((particle) => {
// find closest particle
let closestDist = Infinity
blobs.forEach((blobB) => {
if (this === blobB) return
// optional: ignore far-away blobs here for performance
blobB.particles.forEach((particleB) => {
const dist = particle.distanceTo(particleB)
if (dist < closestDist) closestDist = dist
})
})
// calculate the vector between particle and centre
temp.set(particle.x, particle.y)
temp.sub(this.centre)
temp.normalise()
temp.multiplyScalar(
map(
closestDist,
settings.stoppingDistance[0],
settings.stoppingDistance[1],
0,
0.2,
true
)
)
particle.add(temp)
})
}
draw() {
blobs.forEach((blob) => {
c.beginPath()
c.moveTo(blob.particles[0].x, blob.particles[0].y)
blob.particles.slice(1).forEach((particle) => {
c.lineTo(particle.x, particle.y)
})
c.closePath()
c.fill()
})
}
}
```

When left to run this does just produce a voronoi diagram with a border, but can make more organic shapes if stopped early. By using a variable stopping distance (which ramps down the speed of particles as they approach the minimum distance between blobs) we have finer control over their expansion which makes them look much more natural.

This is starting to feel close to my glass pieces but the blobs still feel quite blocky, like their shapes are entirely derived from their surroundings rather than having any kind of autonomy. Of course they don’t actually have autonomy, they’re code blobs, but they lack a sort of tension which would speak to individuality.

I want the blobs to expand and respond to their surroundings, but also hold their shape in some way. Currently each particle around the perimeter of a blob is acting independently which isn’t realistic if it were a real expanding blob with real perimeter particles. And so we get to soft-body dynamics, a field of computer graphics concerned with making realistic simulations by using real-world physics.

### Mass-spring model

The first thing I added was a spring system; this video by Daniel Schiffman was very instructive on the topic and it’s not very complicated to implement. Essentially every pair of particles on the perimeter is connected by a spring that has a set resting length and tension; if the particles are further away than the spring’s resting length the particles are pulled towards each other, if they’re too close they’re pushed away from one another. The tension describes how much force is applied to get the particles to their resting distance. All the particles in our system have the same mass so we don’t have to take mass into account.

## Code snippet

```
import * as THREE from 'three'
const TENSION = 0.002
const temp = new THREE.Vector2()
class Particle {
constructor(x, y) {
this.pos = new THREE.Vector2(x, y)
this.vel = new THREE.Vector2()
this.acc = new THREE.Vector2()
}
applyForce(force) {
this.acc.plusEq(force)
}
update() {
this.vel.multiplyEq(0.99) // sluggish
this.vel.plusEq(this.acc)
this.pos.plusEq(this.vel)
this.acc.multiplyEq(0)
}
}
class Spring {
constructor(particleA, particleB, restLength) {
this.particleA = particleA
this.particleB = particleB
this.restLength = restLength
}
update() {
let force = temp.copy(this.particleA.pos).sub(this.particleB.pos)
let x = force.length() - this.restLength
force.normalize()
force.multiplyScalar(TENSION * x)
this.particleB.applyForce(force)
force.multiplyScalar(-1)
this.particleA.applyForce(force)
}
}
```

```
import * as THREE from 'three'
const TENSION = 0.002
const temp = new THREE.Vector2()
class Particle {
constructor(x, y) {
this.pos = new THREE.Vector2(x, y)
this.vel = new THREE.Vector2()
this.acc = new THREE.Vector2()
}
applyForce(force) {
this.acc.plusEq(force)
}
update() {
this.vel.multiplyEq(0.99) // sluggish
this.vel.plusEq(this.acc)
this.pos.plusEq(this.vel)
this.acc.multiplyEq(0)
}
}
class Spring {
constructor(particleA, particleB, restLength) {
this.particleA = particleA
this.particleB = particleB
this.restLength = restLength
}
update() {
let force = temp.copy(this.particleA.pos).sub(this.particleB.pos)
let x = force.length() - this.restLength
force.normalize()
force.multiplyScalar(TENSION * x)
this.particleB.applyForce(force)
force.multiplyScalar(-1)
this.particleA.applyForce(force)
}
}
```

This had the effect I was hoping for; particles that would have stretched out into corners are pulled inwards and the whole shape is softened. Over time though the tension causes the blobs to retract to their initial size and eventually fly off the canvas due to the repulsion forces.

### Pressure model

I wanted to add some sort of tension to hold the shape steady once the blob had expanded. Soft-body simulations often use a lattice of springs to achieve this (another Coding Train video covers this comprehensively), but for my blob shapes any springs between the perimeter points would force them back into a circle.

An alternative is to use a pressure-based model. Each blob has a target area, and on each update if the blob’s area is lower than the target area the perimeter particles are pushed out, and if the blob’s area is higher than the target area the perimeter particles are pulled in – it’s the 2D equivalent of a spring.

Implementing this also allowed me to initially increase the size of the blobs by increasing the target area rather than explicitly moving the particles away from the centre of each blob. The added bonus here is that each particle can expand perpendicular to its tangent rather than from the initial seed point.

## Code snippet

```
import * as THREE from 'three'
const temp = new THREE.Vector2()
const origin = new THREE.Vector2()
class Blob {
constructor(x, y, edgePointCount, startSize) {
this.centre = new THREE.Vector2(x, y)
this.maxRadius = startSize
this.edgePointCount = edgePointCount
this.particles = []
this.initialArea = Math.PI * startSize * startSize
this.targetArea = this.initialArea
for (let i = 0; i < this.edgePointCount; i++) {
const angle = (i / this.edgePointCount) * Math.PI * 2
const x = Math.cos(angle) * this.maxRadius + this.centre.x
const y = Math.sin(angle) * this.maxRadius + this.centre.y
this.particles.push(new Particle(x, y))
}
}
update(blobs) {
this.repelBlobs(blobs)
this.maintainPressure()
this.grow()
this.particles.forEach((particle) => {
particle.update()
})
this.springs.forEach((spring) => {
spring.update()
})
}
grow() {
// slowly grow by 3x
if (this.targetArea < this.initialArea * 3) {
this.targetArea *= 1.005
}
}
// https://stackoverflow.com/a/62324338
get area() {
let total = 0
for (let i = 0; i < this._edgePointCount; i++) {
const { x: x1, y: y1 } = this._particles[i]._pos
const { x: x2, y: y2 } =
this._particles[(i + 1) % this._edgePointCount]._pos
total += x1 * y2 - x2 * y1
}
return Math.abs(total / 2)
}
maintainPressure() {
const pressure = this.targetArea / this.area
const forceSize = (pressure - 1) * 20
this.particles.forEach((particle, i) => {
const left =
this.particles[(i + this.edgePointCount - 1) % this.edgePointCount]
const right = this.particles[(i + 1) % this.edgePointCount]
const force = temp
force.copy(left.pos)
force.sub(right.pos)
// perpendicular to particle by averaging position of neighbours
force.rotateAround(origin, Math.PI / 2)
force.normalize()
force.multiplyScalar(forceSize)
particle.applyForce(force)
})
}
}
```

```
import * as THREE from 'three'
const temp = new THREE.Vector2()
const origin = new THREE.Vector2()
class Blob {
constructor(x, y, edgePointCount, startSize) {
this.centre = new THREE.Vector2(x, y)
this.maxRadius = startSize
this.edgePointCount = edgePointCount
this.particles = []
this.initialArea = Math.PI * startSize * startSize
this.targetArea = this.initialArea
for (let i = 0; i < this.edgePointCount; i++) {
const angle = (i / this.edgePointCount) * Math.PI * 2
const x = Math.cos(angle) * this.maxRadius + this.centre.x
const y = Math.sin(angle) * this.maxRadius + this.centre.y
this.particles.push(new Particle(x, y))
}
}
update(blobs) {
this.repelBlobs(blobs)
this.maintainPressure()
this.grow()
this.particles.forEach((particle) => {
particle.update()
})
this.springs.forEach((spring) => {
spring.update()
})
}
grow() {
// slowly grow by 3x
if (this.targetArea < this.initialArea * 3) {
this.targetArea *= 1.005
}
}
// https://stackoverflow.com/a/62324338
get area() {
let total = 0
for (let i = 0; i < this._edgePointCount; i++) {
const { x: x1, y: y1 } = this._particles[i]._pos
const { x: x2, y: y2 } =
this._particles[(i + 1) % this._edgePointCount]._pos
total += x1 * y2 - x2 * y1
}
return Math.abs(total / 2)
}
maintainPressure() {
const pressure = this.targetArea / this.area
const forceSize = (pressure - 1) * 20
this.particles.forEach((particle, i) => {
const left =
this.particles[(i + this.edgePointCount - 1) % this.edgePointCount]
const right = this.particles[(i + 1) % this.edgePointCount]
const force = temp
force.copy(left.pos)
force.sub(right.pos)
// perpendicular to particle by averaging position of neighbours
force.rotateAround(origin, Math.PI / 2)
force.normalize()
force.multiplyScalar(forceSize)
particle.applyForce(force)
})
}
}
```

Between the pressure pushing a blob outwards, the springs contracting, and the repulsion of other blobs, all the forces reach an equilibrium which fixes the shapes in place.

And there we have the final result! It’s not miles away from my initial voronoi explorations, but the irregular soft curves and subtle convex edges lend an organic feel to each piece, and they really do resemble my glass forms quite accurately.

I hope this deep-dive has given you a good insight into the technical journey behind Vitreous. The Truth exhibition is now open so you can peruse all 128 pieces on the exhibition page. In my next article we’ll be taking these 2D forms and turning them into glass with three.js, see you there!