Boid: Rename velocity to direction_vector

[EwoutH]

A velocity has a magnitude and a direction. Since the velocity was modelled as a unit vector, the magnitude is constant, and it only represents the direction, not the velocity (which is magnitude + direction). Therefor, rename velocity to direction_vector.

In this model, the magnitude is represented by speed.

Also removed a redundant division.

[rht]

The / 2 affects the result. This is self.velocity +=, not self.velocity =.

[EwoutH]

Good catch, you're totally right. I will reply in the main thread more extensively.

[tpike3]

Ohhh.. just second I was reading this background... whether the / 2 is need is the debate correct?

[rht]

I agree with the rename, but the removal of the division by 2 has a problem. To rephrase what I said, adding a unit vector with another vector with a different magnitude (with or without division by 2) and then to normalize it to unit vector again will result with a unit vector with a different direction.

[EwoutH]

You're right on the updating and not overwriting, I read over that.

However:

• is updating the desired behavior? Should they update their direction every step (like now) or should they just create a new direction anyways?
• Where does the / 2 come from? Is it correcting for that there are three factors (shouldn't it be / 3). Or is it some weight? In that case, * 0.5 is more logical, and it should be maybe an input variable (like update_direction_weight).

The / 2 feels like a magic number. Do we know where it came from and why?

[rht]

is updating the desired behavior?

Updating makes more sense, because it implements an acceleration. There should be an inertia in the boid movement.

Where does the / 2 come from?

Not sure. https://www.red3d.com/cwr/boids/ has the original implementation

(defmethod (:NAVIGATE BOID-BRAIN) ()
  (vlet* ((avoid-obstacles   (send self :environmental-obstacle-avoidance))
	  (avoid-flockmates  (send self :flockmate-avoidance))
	  (velocity-matching (send self :velocity-matching))
	  (centering-urge    (send self :neighborhood-centering))
	  (migratory-urge    (send self :migratory-urge))
	  (course-leveling   (send self :course-leveling))
	  (course-damping    (send self :course-damping)))
    ;;
    ;; Available-acceleration should probably be 1.0, but I've set it a little higher to
    ;; avoid have to readjust the weighting factors for all of the acceleration requests.
    ;;
    (vlet* ((composite-acceleration (prioritized-acceleration-allocation 1.3 ;; 1.0
									 avoid-obstacles
									 avoid-flockmates
									 velocity-matching
									 centering-urge
									 migratory-urge
									 course-leveling
									 course-damping))
	    (weights                (values-list (send self :steering-weights)))
	    (desired-acceleration   (componentwise-product weights composite-acceleration)))
      '(when (eql self (default-boid))
	 (format w " n"))
      (values-list desired-acceleration))))

(defun PRIORITIZED-ACCELERATION-ALLOCATION (available-acceleration &rest vectors)
  (loop with allocated-acceleration = 0
	with (tx ty tz) = '(0 0 0)
	for v in vectors
	for (vx vy vz) = (or v '(0 0 0))
	for mag-v = (magnitude v)
	for scaler = (if (> (+ allocated-acceleration mag-v) available-acceleration)
			 (// (- available-acceleration allocated-acceleration) mag-v)
			 1)
	do (incf tx (* vx scaler))
	   (incf ty (* vy scaler))
	   (incf tz (* vz scaler))
	   (incf allocated-acceleration mag-v)
	   ;; (format tow "~&  (~d ~D ~d), ~d, ~d" tx ty tz allocated-acceleration mag-v)
	until (> allocated-acceleration available-acceleration)
	finally
	  ;; (format t "~&  mag= ~d" (magnitude-xyz tx ty tz))
	  ;; (format w "~&~D" allocated-acceleration)
	  (return (values tx ty tz))))

update_direction_weight makes sense to me. But it seems to differ from available-acceleration in that code.

[EwoutH]

Interesting. If we’re actually talking about acceleration, I find it incredibly weird that it get’s normalized to a unit vector each step.

Can we find the original model conceptualization or formalization?

[rht]

Superseded by #104, except for the discussion.

@coderbeta1 if you want to read the original implementation of Boid Flockers, the thread above has a link and relevant code snippet to it. We might have deviated from the original implementation by not implementing acceleration.

[coderbeta1]

Proposed Implementation:

This seems like a good implementation of the boid flocker model with acceleration implemented as intended.

TLDR:

1. Boids: The boids have a self.position vector, self.velocity vector and self.acceleration vector.
2. Position Update: new self.position += self.velocity
3. Velocity Update: new self.velocity += self.acceleration
4. Max Speed: To prevent too high magnitude of self.velocity, there is a limiter. Method for limiter is
   if || self.velocity || > max_speed , then get unit vector in direction of velocity and multiple by max_speed. (So now max speed of boid is max_speed)
5. Acceleration Update: There are 3 things that affect the acceleration

1. Alignment Force: Steer towards average direction in neighbourhood. How?

• alignment = avg_velocity (normalize and multiply with max_speed) -self.velocity

2. Separation Force: Move boids away so they don't crash. How?

• Find average_dist vector, multiply by -1 (to move in opposite direction).
• repulsion = average_dist (normalize and multiply with max_speed) - self.velocity.
• There is a max_force.
• if || repulsion || > max_force , then get unit vector in direction of repulsion and multiple by max_force. (So now max steering of boid is max_force)

3. Cohesion Force: Bring the boids together. How?

• vec_to_com = center_of_mass - self.position
• cohesion = vec_to_com (normalize and multiply with max_speed) - self.velocity.
• There is a max_force.
• If || cohesion || > max_force , then get unit vector in direction of steering and multiple by max_force. (So now max steering of boid is max_force)

New self.acceleration += (alignment + cohesion + repulsion)

Current Implementation:

• Concept of Acceleration is not implemented: This means that normal forces can't be applied. The boid will directly try to move towards final acceleration. This makes the implementation look not so smooth
• Concept of Centre of Mass is not Implemented: Cohesion Force is implemented incorrectly in my opinion. It just takes vector between two boids and adds them. Reason - I don't know
• Separation Force: It seems fine, but has a minimum separation factor - Why? I don't know
• Match-Vector/Alignment Force: Seems fine. But uses net direction alone, i.e. a velocity is applied instead of acceleration

[EwoutH]

Thanks for working on this!

I think the most important distinction that needs to be made is if you want to use a stationary coordinate system, of if you want to use the boid itself as the frame of reference.

Both have advantages. I'm nudging towards the later, because there is a good case to be made that boids have more control for changing their direction (heading) than for changing their velocity (speed). For me it also feels more intuitive.

[coderbeta1]

a good case to be made that boids have more control for changing their direction (heading) than for changing their velocity (speed)

True that. This just means that with the proposed implementation we increase the max_force. This would mean that the boids are more capable of changing direction (since top speed is limited).

If we talk about intuition I think it completely depends on the person.

For a third person viewer the stationary coordinate system is more intuitive. This is because

• it is easy to see the cluster and imagine things that way.
• Concept of Centre of Mass (COM) is more intuitive when we see the whole picture.

On the other hand,

• Using the frame of reference of the boid will let us apply changes easily but it is a little more difficult to imagine change in direction towards COM.
• Using an accelerating frame of reference and imaging a pseudo force is more convenient

Calculating the forces is easy to implement either way. I feel that the way the article is implemented is quite intuitive when is comes to acceleration.

[rht]

I agree with proceeding with the implementation with acceleration and center of mass. Let's proceed forward. The citation (pseudocode at https://vergenet.net/~conrad/boids/pseudocode.html) at the README.md also seems to suggest a form of acceleration.

This is the original PR that implemented the example: projectmesa/mesa#136, which has no information about the source of the model.

[coderbeta1]

ok so I tried to follow the pseudocode as mentioned by @rht . After spending a lot of time reimplementing it, what I noticed is that the original implementation is the same as https://vergenet.net/~conrad/boids/pseudocode.html . Its just that the code implements it in a way that is not easy to understand (since there are no comments in the logic part).

• The part of centre of mass is actually implemented very cleverly already. (It is just implemented in a different frame of reference - with boid as the origin)
• In the citation they implement static values of scaling factor for cohesion, separation and alignment. But with the current implementation it is easy to fine tune and check the importance of the factors

One thing that I feel should be changed is the scheduling.

• For this case specifically I think Simultaneous Activation is more logical than random. I tried to implement it and the result seems to be better with simultaneous activation. It just feels more smooth and predictable.
• If we consider a cluster of boids with random activation and a cluster with simultaneous, it just seems logical to go with simultaneous. This is because if we take an example of (lets say) 3 boids and do random activation, the order of steps will make a huge difference (and in some cases compound the result)

[rht]

@coderbeta1 you are doing a scholarly work here! We should update Agents.jl about your findings once the dust has settled, as they use Boid for https://github.com/JuliaDynamics/ABMFrameworksComparison.

I have a question regarding with the acceleration: the v1, v2, v3, ... in https://vergenet.net/~conrad/boids/pseudocode.html are technically the forces / mass * deltat. But in the current mesa-examples implementation, the speed is held constant, and only the direction changes. Does the acceleration part differ from the citation?

I agree with the SimultaneousActivation being more logical. Eliminating first mover advantage and order hysteresis makes sense.

[EwoutH]

Awesome work so far!

For this case specifically I think Simultaneous Activation is more logical than random.

Agreed! It's also faster and less noisy, and all agents act in continuous time independent from each other, so it's obviously the way to go. RandomActivation also allows an agent acting two times in a row, which from a reaction-time standpoint is not the proper way to model it.

[coderbeta1]

I have a question regarding with the acceleration: the v1, v2, v3, ... in https://vergenet.net/~conrad/boids/pseudocode.html are technically the forces / mass * deltat. But in the current mesa-examples implementation, the speed is held constant, and only the direction changes. Does the acceleration part differ from the citation?

@rht So the thing is, the current implementation does differ from the citation. It does not have acceleration as such:

• On implementing it exactly as the paper: Issue is that the visualization just does not converge. This is due to exploding/vanishing velocity and acceleration.
• On deviating a bit and implementing max speed and acceleration: This looks amazing and is almost perfect, but the issue is that it again does not converge - the reason being that there are way too many hyperparameters to fine tune, namely cohere_factor, align_factor, separation_factor, max_speed, max_force and start_speed. (I have some values which somewhat work. I just chose them by trial and error)

So while implementing acceleration with clipping max speed and acceleration is way better than without, it is very difficult to figure out what values to choose to get best representation of real life scenario.

With the current implementation the solution does at the very least, converge

If the aim is to leave it to the user to choose the values and understand the factors which affect the simulation, I think we should go with acceleration and max_speed and max_force

[rht]

If you look at Reynold's original implementation in Lisp, https://www.red3d.com/cwr/code/boids.lisp, it has min-speed, max-speed, and max-accel-relative. But it also has acceleration-smoothing, a gravitational field, etc. It's too detailed, but maybe we could have a folder boid_flockers_full_version later on.

If you look at the pseudocode in https://vergenet.net/~conrad/boids/pseudocode.html, it does have a speed limit. Maybe just limiting the speed would be sufficient, without having to set an upper bound to the acceleration?

[coderbeta1]

Apologies for the late reply, got some important work on my side.

If you look at the pseudocode in https://vergenet.net/~conrad/boids/pseudocode.html, it does have a speed limit. Maybe just limiting the speed would be sufficient, without having to set an upper bound to the acceleration?

The speed limit does help, but it does not work if acceleration does not have a limit.

• In general I noticed that with max acceleration < max speed I got best results
• Changing the vision radius of boid has a massive difference on the solution. High radius gives best results (it does expectedly slow down the simulation though)

It's too detailed, but maybe we could have a folder boid_flockers_full_version later on.

• Seems like a great idea.
• I think for now I'll just create a PR for random sequence to simultaneous and the max_speed + max_acceleration change