fix: lerp target in ParamChange::MLinAlg

Pressing 'c' will center the screen around the average position of the
particles. Pressing 'f' will center the screen on each frame.

Pressing 'b' will center the screen around the best found point.

.gitignore

@@ -1,3 +1,4 @@
 .clangd
 *.o
 bin/swarm
+obj/

Makefile

@@ -6,19 +6,16 @@ OBJ_DIR = ./obj
 
 HEADERS   = $(wildcard $(INC_DIR)/*.hpp)
 SRC_FILES = $(wildcard $(SRC_DIR)/*.cpp)
+DEP_FILES = $(OBJ_FILES:.o=.d)
 OBJ_FILES = $(addprefix $(OBJ_DIR)/,$(notdir $(SRC_FILES:.cpp=.o)))
 
-CXXFLAGS += -Wall -std=c++20 -O2
+CXXFLAGS += -Wall -std=c++20 -O2 -I $(INC_DIR)
 
 all: bin/swarm
 
-$(OBJ_DIR)/%.o: src/%.cpp
-	$(CXX) -I $(INC_DIR) $(CXXFLAGS) -c -o $@ $^
-
 bin/swarm: $(OBJ_FILES)
-	@mkdir -p $(OBJ_DIR)
 	@mkdir -p bin
-	$(CXX) -I $(INC_DIR) $(CXXFLAGS) -o $@ $^
+	$(CXX) $(CXXFLAGS) -o $@ $^
 
 .PHONY: debug
 debug: CXXFLAGS += -O0 -g
@@ -28,3 +25,9 @@ debug: clean bin/swarm
 clean:
 	$(RM) obj/*.o
 	$(RM) bin/swarm
+
+-include $(DEP_FILES)
+
+$(OBJ_DIR)/%.o : src/%.cpp
+	mkdir -p $(@D)
+	$(CXX) $(CXXFLAGS) -MMD -c $< -o $@

README.md

@@ -0,0 +1,15 @@
+# swarmc
+
+A [Particle Swarm
+Optimisation](https://en.wikipedia.org/wiki/Particle_swarm_optimization)
+visualiser in C++. Run with `swarm <iteration count>` and enjoy the pretty
+colors. Pressing `h` will print out a help menu.
+
+Requires an ANSI terminal (emulator) with support for 24-bit colours to render.
+
+The implemented showcase works for `R^2 -> R` functions, but the underlying code
+can work for much more, including for spaces other than real spaces (see the
+type parameters in `Agent` et al).
+
+Supports dynamically changing Particle parameters.
+Although not implemented (yet!), the code supports arbitrary swarm topologies

include/screen.hpp

@@ -40,19 +40,23 @@ struct Screen {
     Symbol &at(int x, int y);
     Symbol &at(float x, float y);
 
-    void resize(std::size_t n);
     void resize(std::size_t x, std::size_t y);
 
+    void move_to(float x, float y);
+
     std::pair<float, float> screen_to_xy(int w, int h);
     std::pair<int, int> xy_to_screen(float x, float y);
 
-    Screen(std::size_t n) : buf(nullptr) { resize(n); }
-    Screen(std::size_t x, std::size_t y) : buf(nullptr) { resize(x, y); }
+    Screen(std::size_t x, std::size_t y) : buf(nullptr)
+    { resize(x, y); move_to(0, 0); }
     ~Screen() { delete[] buf; }
 
     std::vector<const vec<2>*> points;
+    std::vector<std::pair<const vec<2>*,const vec<2>*>> vecs;
     char_shader_t shader;
 
+    bool draw_vecs = true;
+
 private:
     Symbol *buf;
     Symbol _dummy;

include/swarm.hpp

@@ -1,5 +1,6 @@
 #pragma once
 
+#include <cmath>
 #include <cstddef>
 
 #include <functional>
@@ -7,13 +8,37 @@
 #include <vec.hpp>
 #include <vector>
 
+/**
+ * An Agent in a Particle Swarm Optimisation simulation.
+ * #T is the type of the position vector the Agent uses.
+ */
 template<typename T>
 struct Agent {
+    /**
+     * Type of the fitness function.
+     */
     using F = std::function<float(T)>;
 
+    /**
+     * Moves the Agent along its velocity by a factor of #dt.
+     *
+     * The parameter #dt is useful when drawing the Agent moving faster than it
+     * #tick()s.
+     *
+     * @see: step()
+     */
     virtual void move(float dt) = 0;
+
+    /**
+     * Steps the algorithm values one tick forward.
+     * Warning: Does not move the Agent. Use #move() for this
+     * @see: move()
+     */
     virtual void step() = 0;
 
+    /**
+     * The best position this Agent has seen during its lifetime.
+     */
     virtual std::pair<float, T> best() const = 0;
 
     Agent(F f) : f(f) {}
@@ -23,6 +48,8 @@ protected:
 
 template <typename A>
 concept ParameterChangeAlgorithm = requires(A alg, float f, unsigned int i) {
+    // initial viscosity, initial nostalgia, initial peer pressure, current
+    // iteration, max iteration
     { alg(f, f, f, i, i) } -> std::same_as<vec<3>>;
 };
 
@@ -39,9 +66,23 @@ namespace ParamChange {
 	    float visc, float nostal, float peerp, int curr, int max) {
 	    float pct_done = (float)curr / (float)max;
 	    return {
-		visc   - (0.9f - 0.4f) * pct_done,
-		nostal - (2.5f - 0.5f) * pct_done,
-		peerp  + (2.5f - 0.5f) * pct_done,
+		visc   - (visc   - 0.4f) * pct_done,
+		nostal - (nostal - 0.5f) * pct_done,
+		peerp  - (peerp  - 2.5f) * pct_done,
+	    };
+	}
+    };
+
+    template<auto steepness = 3>
+    struct MTanhAlg {
+	vec<3> operator()(
+	    float visc, float nostal, float peerp, int curr, int max) {
+	    float pct_done = (float)curr / (float)max;
+	    float scale_factor = (std::tanh(steepness * (pct_done - 0.5)) + 1) / 2;
+	    return {
+		visc   - (visc   - 0.4f) * scale_factor,
+		nostal - (nostal - 0.5f) * scale_factor,
+		peerp  - (peerp  - 2.5f) * scale_factor,
 	    };
 	}
     };
@@ -54,6 +95,16 @@ struct Particle : public Agent<vec<N>> {
     void move(float dt=1) override { position = position + velocity * dt; }
 
     void step() override {
+	// Get the new parameter values
+	vec<3> params = alg(kviscosity, knostalgia, kpeer_pressure,
+			    curr_iter, max_iter);
+
+	viscosity = params.x;
+	nostalgia = params.y;
+	peer_pressure = params.z;
+
+	curr_iter++;
+
 	velocity = viscosity * velocity
 	    + nostalgia * (pb_pos - position)
 	    + peer_pressure * (peer.best().second - position);
@@ -63,15 +114,6 @@ struct Particle : public Agent<vec<N>> {
 	    pb = y;
 	    pb_pos = position;
 	}
-
-	vec<3> params = alg(kviscosity, knostalgia, kpeer_pressure,
-			    curr_iter, max_iter);
-
-	viscosity = params.x;
-	nostalgia = params.y;
-	peer_pressure = params.z;
-
-	curr_iter++;
     };
 
     const float kviscosity = 0.9f;
@@ -90,6 +132,7 @@ struct Particle : public Agent<vec<N>> {
     }
 
     const vec<N> &get_position() const { return position; };
+    const vec<N> &get_velocity() const { return velocity; };
 
 private:
     A alg;
@@ -146,7 +189,7 @@ struct Swarm : public Agent<vec<N>> {
 
     const std::vector<Particle<N>>& get_particles() { return particles; };
 
-    Swarm(Agent<vec<N>>::F f, unsigned max_iter = 150)
+    Swarm(Agent<vec<N>>::F f, unsigned max_iter)
 	: Agent<vec<N>>(f), max_iter(max_iter) {}
 private:
     unsigned max_iter;

src/main.cpp

@@ -1,4 +1,4 @@
-#include "screen.hpp"
+#include <screen.hpp>
 #include "vec.hpp"
 #include <cmath>
 #include <cstdio>
@@ -8,6 +8,7 @@ static constexpr int   kFPS = 60;
 static constexpr float kDT  = 1.0 / kFPS;
 
 float f(vec<2> x) {
+    // return std::pow(x.x, 2) + std::pow(x.y, 2);
     return  50 * (std::pow(std::sin((x.x - 10)/2),2) + std::pow(std::sin(x.y/2),2)) +
 	std::pow(std::abs(x.x - 3.141592653589), 1.2) + std::pow(std::abs(x.y), 1.2);
 }
@@ -16,7 +17,13 @@ int
 main(int argc, char **argv) {
     int iter_count;
     if(argc < 2) iter_count = 100;
-    else sscanf(argv[1], "%d", &iter_count);
+    else {
+	int scanned = sscanf(argv[1], "%d", &iter_count);
+	if(!scanned) {
+	    printf("Usage: %s <iteration count>\n", argv[0]);
+	    exit(1);
+	}
+    }
 
     Swarm<2> swarm(f, iter_count);
     Screen scr(80, 24);
@@ -47,6 +54,10 @@ main(int argc, char **argv) {
 
     for(const auto &p : swarm.get_particles()) {
 	scr.points.push_back(&p.get_position());
+	scr.vecs.push_back({
+		&p.get_position(),
+		&p.get_velocity(),
+	    });
     }
 
     enter_noncanonical_mode();
@@ -55,6 +66,7 @@ main(int argc, char **argv) {
 
     bool pause = false;
     bool frame_step = false;
+    bool auto_follow = false;
 
     // initialize colorizer scale/translate
     scr.draw();
@@ -63,21 +75,46 @@ main(int argc, char **argv) {
 
     auto update_and_draw = [&]() -> void {
 	scr.draw();
+	auto [y, x] = swarm.best();
+	printf("\033[30;40m\33[2K\r\033[48;2;%d;%d;%dm"
+	       "Current best: (%f, %f) = %f\n",
+	       255, 255, 255, x.x, x.y, y);
     };
 
+    auto center_screen = [&]() -> void {
+	vec<2> pos = 0;
+	for(const auto &p : swarm.get_particles()) {
+	    pos += p.get_position();
+	}
+	pos /= swarm.get_particles().size();
+	scr.move_to(pos.x, pos.y);
+    };
+
+    // We draw to the screen at a rate of `kFPS', but step()ing the swarm at
+    // this rate would be far too fast to be interesting to look at. On the
+    // other hand, step()ing once a second is too slow.
+    //
+    // So, we step() every quarter second. To line up the move() so that a unit
+    // move() happens preceding each step(), we also need to move as fast as
+    // we've shortened our step() interval.
     for(int i = 0; i < iter_count * kFPS/4;) {
 	if(!pause) {
 	    scr.clear();
 
 	    update_and_draw();
 
+	    // Since we step() 4x a second, we need to move() 4x as fast as moving by dt
+	    // each frame.
 	    swarm.move(kDT * 4);
 	    if(i % (kFPS/4) == (kFPS/4)-1)
 		swarm.step();
 
 	    ++i;
+	    printf("Current iteration: %d\nCurrent frame: %d\n", i * 4 / kFPS, i);
 	}
 
+	if(auto_follow) center_screen();
+
 	if(frame_step) {
 	    pause = true;
 	    frame_step = false;
@@ -91,6 +128,14 @@ main(int argc, char **argv) {
 	case ' ':
 	    pause = !pause;
 	    break;
+	case ',':
+	    for(int j = 0; j < kFPS/4; ++j) {
+		swarm.move(kDT * 4);
+		if(i % (kFPS/4) == (kFPS/4)-1)
+		    swarm.step();
+
+		++i;
+	    }
 	case '.':
 	    pause = false;
 	    frame_step = true;
@@ -146,12 +191,34 @@ main(int argc, char **argv) {
 	    update_and_draw();
 	    break;
 
-	case 'p':
-	    auto [y, x] = swarm.best();
-	    printf("Current best: f(%.10f, %.10f) = %.10f", x.x, x.y, y);
+	case 'v':
+	    scr.draw_vecs = !scr.draw_vecs;
+	    break;
+
+	case 'b': {
+	    const auto [_, b] = swarm.best();
+	    scr.move_to(b.x, b.y);
 	    break;
 	}
 
+	case 'c':
+	    center_screen();
+	    break;
+
+	case 'f':
+	    auto_follow = !auto_follow;
+	    break;
+
+	case 'h':
+	    printf("  movement  zoom  coloring  pause step \n"
+		   "     W       io    IK LO     SPC   .,  \n"
+		   "    ASD                                \n"
+		   "   [%c] draw (v)elocities   quit:  q   \n",
+		   scr.draw_vecs ? 'x' : ' ');
+	    break;
+	}
+
+	// Wait for kDT (time between two frames)
 	usleep(1000 * 1000 / kFPS);
     }
 
@@ -161,7 +228,7 @@ cleanup:
     enter_canonical_mode();
 
     auto [y, x] = swarm.best();
-    printf("Best: f(" V2_FMT ") = %.3f", V2_ARG(x), y);
+    printf("Minimum found:\n f(%f, %f) = %f\n", x.x, x.y, y);
 
     return 0;
 }

src/screen.cpp

@@ -4,11 +4,6 @@
 
 /* Screen */
 
-void Screen::resize(std::size_t n) {
-    w = n; h = 1;
-    buf = (Symbol *)realloc(buf, n * sizeof(Symbol));
-}
-
 void Screen::resize(std::size_t x, std::size_t y) {
     w = x; h = y;
     buf = (Symbol *)realloc(buf, x * y * sizeof(Symbol));
@@ -41,6 +36,10 @@ std::pair<int, int> Screen::xy_to_screen(float x, float y) {
     };
 }
 
+void Screen::move_to(float x, float y) {
+    dx = x - w*sx/2;
+    dy = y - h*sy/2;
+}
 
 void Screen::clear() {
     static const Symbol s {' ', {0,0,0}};
@@ -57,6 +56,32 @@ void Screen::draw() {
 	}
     }
 
+    // draw lines for each velocity vector
+    if(draw_vecs) {
+	for(const auto &[start, vel] : vecs) {
+	    auto end = *start + (*vel) / 3;
+	    if(start->x == end.x) {
+		// TODO: draw horizontal/vertical lines
+		continue;
+	    }
+
+	    auto a = (end.y - start->y) / (end.x - start->x);
+	    auto b = start->y - a * start->x;
+
+	    for(int i = 0; i < w; ++i) {
+		const auto [x, _] = screen_to_xy(i, 0);
+		if(x < std::min(start->x, end.x)
+		   || x > std::max(start->x, end.x)) continue;
+
+		auto y = a*x + b;
+		at(x, y) = Symbol{
+		    .sym = '+',
+		    .color = { 0.75, 0, 0 },
+		};
+	    }
+	}
+    }
+
     // write out a '#' wherever we have a point registered
     for(const auto &p : points) {
 	at(p->x, p->y) = {.sym = '#', .color = 1};