Robotic swarm simulation: aggregation

Question

This program simulates the phenomenon of aggregation of a swarm of robots with no cooperation between each other and each robot following a very simple rule:

The velocity of this robot has a random direction and the absolute value of this velocity is inversely proportional to the number of robots inside of a threshold distance from this robot.

This means that the robot will move slower the more friends it has nearby, so it will spend more time near other robots rather than far away, and this is enough to create aggregation of robots, with no communication nor cooperation needed.

The program does a live representation of the swarm of robots and a live plot of the average number of neighbors of each robot, that is increasing over time.

Screenshot after some time

Graph of average number of neighbors over time

Graph of average number of neighbors over time, higher time limit

Code

import sys
import pygame
from pygame.locals import QUIT
import random
import math
from operator import mul, truediv
from statistics import mean
from itertools import combinations
import matplotlib.pyplot as plt
import numpy as np


class Robot:
    def __init__(self, x, y, v=0, v0=None, aggregating_mode=True):
        self.x = x
        self.y = y
        self.v = 0
        self.aggregating_mode = aggregating_mode
        self.v0 = v0 or (20 if aggregating_mode else 2)

    def distance_to_robot(self, other):
        return math.sqrt((self.x - other.x) ** 2 + (self.y - other.y) ** 2)

    def count_nears(self, robots, threshold=100):
        return sum(1 for robot in robots if self.distance_to_robot(robot) < threshold)

    def update_speed(self, robots, threshold=100):
        # Denominator is never 0 because the robot will always count at least itself
        operation = truediv if self.aggregating_mode else mul
        self.v = self.v0 * operation(1.0, self.count_nears(robots, threshold=threshold))

    def update_position(self, width, height):
        theta = random.random() * 2 * math.pi
        self.x += math.cos(theta) * self.v
        self.y += math.sin(theta) * self.v
        self.x %= width
        self.y %= height

    def draw_self(self, SCREEN):
        color = (255, 0, 0) if self.aggregating_mode else (0, 255, 0)
        pygame.draw.circle(SCREEN, color, (self.x, self.y), 10)

    def __repr__(self):
        return f"Robot({self.x}, {self.y}, {self.v}"


def main():
    pygame.init()
    pygame.display.set_caption("Swarm robotics simulation: aggregation and distancing")
    FPS = 60
    FPS_CLOCK = pygame.time.Clock()

    pygame.display.set_mode((0, 0), pygame.FULLSCREEN)
    WIDTH, HEIGHT = pygame.display.Info().current_w, pygame.display.Info().current_h
    SCREEN = pygame.display.set_mode((WIDTH, HEIGHT))

    def random_pos():
        return (random.randint(0, WIDTH), random.randint(0, HEIGHT))

    robots = [Robot(*random_pos(), aggregating_mode=True) for _ in range(140)]
    # Game loop.
    i = 0
    average_counts = []
    is_ = []
    PLOTTING = True
    while True:
        i += 1
        SCREEN.fill((0, 0, 0))

        for event in pygame.event.get():
            if event.type == QUIT:
                pygame.quit()
                sys.exit()

        # Update and Draw
        for robot in robots:
            robot.update_speed(robots)
            robot.update_position(WIDTH, HEIGHT)
            robot.draw_self(SCREEN)

        if i % 100 == 0:
            average_near_count = mean((robot.count_nears(robots) for robot in robots))
            average_counts.append(average_near_count)

            if PLOTTING:
                is_.append(i)
                plt.axis([0, 15000, 0, 14])
                plt.xlabel("Timestep")
                plt.ylabel("Average number of neightbours")

                plt.scatter(is_, average_counts)
                plt.pause(0.1)

                plt.show()

            print(average_near_count)
        pygame.display.flip()
        FPS_CLOCK.tick(FPS)


if __name__ == "__main__":
    main()

Answer 1

Your user interface is fairly awkward, halting the simulation periodically to show and then destroy a graph.

You've chosen to use pygame, but matplotlib has built-in animation support. I am going to suggest that you switch to that, because:

• It will allow for a seamless, concurrent display of the bot positions and whatever other statistics you want, in separate portions of a figure grid
• This is an analysis project and not strictly a game

Your representation of a single robot as a class, whereas it's mathematically sound, is inefficient. Numpy is good at vectorisation, which means

• Don't represent a robot as a single object; represent it as a row in an array
• Don't even separate your x and y coordinates; keep them in columns of one array
• Don't call count_nears once for every robot; instead set up a three-dimensional matrix: source robots by destination robots by x/y, and apply the Euclidean norm over the outermost axis to get a vectorised distance.
• Don't call into math or operator; stick to Numpy.

Spelling: neightbours -> neighbours.

I find a scatterplot for your mean to be an odd choice, and would expect a line plot instead.

Topology

I'll have to ask for forgiveness from my mathematician friends, since I'm speaking out of my depth:

You've (perhaps inadvertently) created a coordinate space topologically homeomorphic to a toroid. This is due to your modulus wrap-around. This is not bad in itself, but it's incompatible with another decision of yours - you have a discontinuity in your use of the Euclidean norm.

Imagine two points near the left-hand side. Their distance is low. Then one point drifts past the boundary and moves to the far right-hand side. Their distance is high. But should it be, really? The amount of space one point needs to traverse to intersect with the other point is still low, so the distance should be considered low.

There are two different ways to address this. To keep your space fully continuous, you'll need to change your use of the Euclidean norm such that - for each point - you don't measure against one other point; you measure against the shortest toroidal path among four.

The less continuous - but certainly the simpler - to implement, is to replace your modulus with max and min calls. I've not implemented either in my suggested code.

Suggested

from typing import Optional

import matplotlib.pyplot as plt
import numpy as np
from matplotlib.animation import FuncAnimation
from matplotlib.collections import Collection
from numpy.random import default_rng


FPS = 60


class Swarm:
    def __init__(
        self,
        width: float = 1000,
        height: float = 1000,
        near_threshold: float = 100,
        v0: Optional[float] = None,
        n: int = 140,
        aggregate: bool = True,
        seed: int = 0,
    ) -> None:
        self.width, self.height, self.near_threshold, self.n, self.aggregate = (
            width, height, near_threshold, n, aggregate,
        )
        self.v0 = v0 or (20 if aggregate else 2)
        self.rand = default_rng(seed=seed).random
        self.pos = width*self.rand((2, n))
        self.n_near_mean: float
        self.size = np.array((width, height))[:, np.newaxis]

    def get_v(self) -> tuple[
       np.ndarray,  # new speeds
       float,       # new mean proximate-neighbour count
    ]:
        displacement = self.pos[:, :, np.newaxis] - self.pos[:, np.newaxis, :]
        distance = np.linalg.norm(displacement, axis=0)
        n_near = np.count_nonzero(distance < self.near_threshold, axis=0)
        if self.aggregate:
            v = self.v0 / n_near
        else:
            v = self.v0 * n_near
        return v, n_near.mean()

    def update(self) -> None:
        v, self.n_near_mean = self.get_v()
        theta = 2*np.pi*self.rand(self.n)
        projections = np.stack((
            np.cos(theta), np.sin(theta),
        ))
        self.pos = np.mod(self.pos + v*projections, self.size)


def display(swarm: Swarm) -> tuple[plt.Figure, FuncAnimation]:
    fig: plt.Figure
    ax_bot: plt.Axes
    ax_stats: plt.Axes
    fig, (ax_bot, ax_stats) = plt.subplots(nrows=1, ncols=2)

    bot_scatter: Collection = ax_bot.scatter(*swarm.pos)
    ax_bot.set_title('Robot positions')
    ax_bot.set_xlabel('x')
    ax_bot.set_ylabel('y')
    ax_bot.axis((0, swarm.width, 0, swarm.height))

    frames = []
    means = []
    stat_line: plt.Line2D
    stat_line, = ax_stats.plot(frames, means)
    ax_stats.set_title('Mean proximate neighbours')
    ax_stats.set_xlabel('Timestep')
    ax_stats.set_ylabel('Count')
    ax_stats.axis((0, 15_000, 0, 14))

    def init() -> tuple[plt.Artist, ...]:
        return bot_scatter, stat_line

    def update(frame: int) -> tuple[plt.Artist, ...]:
        swarm.update()
        bot_scatter.set_offsets(swarm.pos.T)
        if frame % FPS:
            return bot_scatter,

        print(f'{swarm.n_near_mean:.1f}')
        frames.append(frame)
        means.append(swarm.n_near_mean)
        stat_line.set_data(frames, means)
        return bot_scatter, stat_line

    anim = FuncAnimation(
        fig=fig, func=update, init_func=init, interval=1e3/FPS, blit=True,
    )
    return fig, anim


def main() -> None:
    swarm = Swarm()
    print('Mean proximate neighbour count:')
    fig, anim = display(swarm)
    plt.show()


if __name__ == "__main__":
    main()

Output