Cleanup of Optimization file structure

Partition of swarm::mpi_run() into multiple functions for clarity

src/Optimization/CMakeLists.txt

@@ -3,7 +3,12 @@ project(Oersted_Optimization)
 set(SOURCE_FILES
         ./include/Optimization.hpp
 
-        include/Optimization.hpp src/Particle.cpp src/Particle.h)
+        ./src/Particle.cpp ./src/Particle.h
+        ./src/ContiguousDataMap.cpp ./src/ContiguousDataMap.h
+        ./src/CoordinateState.cpp ./src/CoordinateState.h
+        ./src/Coordinate.cpp ./src/Coordinate.h
+        ./src/CoordinateSpace.cpp ./src/CoordinateSpace.h
+        ./src/Swarm.cpp ./src/Swarm.h)
 
 add_library(optimization SHARED ${SOURCE_FILES})
 

src/Optimization/include/Optimization.hpp

 #ifndef OERSTED_OPTIMIZATION_HPP
 #define OERSTED_OPTIMIZATION_HPP
 
+#include "../src/ContiguousDataMap.h"
+#include "../src/Coordinate.h"
+#include "../src/CoordinateSpace.h"
+#include "../src/CoordinateState.h"
 #include "../src/Particle.h"
+#include "../src/Swarm.h"
 
 #endif //OERSTED_OPTIMIZATION_HPP

src/Optimization/src/ContiguousDataMap.cpp

+#include "ContiguousDataMap.h"

src/Optimization/src/ContiguousDataMap.h

+#ifndef OERSTED_CONTIGUOUSDATAMAP_H
+#define OERSTED_CONTIGUOUSDATAMAP_H
+
+#include <vector>
+#include <map>
+#include <cmath>
+
+template <class Key, class T>
+class ContiguousDataMap {
+public:
+    T &operator[] (Key k) {
+        auto it = Map.find(k);
+        if (it != Map.end()) {
+            return Vector[it->second];
+        } else {
+            Map[k] = Vector.size();
+            Vector.push_back(T());
+            return Vector.back();
+        }
+    };
+
+    T &operator[] (size_t k) { return Vector[k]; };
+
+    std::map<Key,size_t> &map() { return Map; };
+
+    std::vector<T> &vector() { return Vector; };
+
+    size_t size() const { return Vector.size(); };
+
+protected:
+    std::map<Key,size_t> Map;
+
+    std::vector<T> Vector;
+};
+
+using ObjectiveMap = ContiguousDataMap<std::string,double_t>;
+
+#endif //OERSTED_CONTIGUOUSDATAMAP_H

src/Optimization/src/Coordinate.cpp

+#include "Coordinate.h"

src/Optimization/src/Coordinate.h

+#ifndef OERSTED_COORDINATE_H
+#define OERSTED_COORDINATE_H
+
+#include <cmath>
+
+class Coordinate {
+public:
+    Coordinate() {};
+
+    Coordinate(double_t lb, double_t nom, double_t ub) : LowerBound{lb}, Nominal{nom}, UpperBound{ub} {};
+
+    double_t lower_bound() const { return LowerBound; };
+
+    double_t nominal() const { return Nominal; };
+
+    double_t upper_bound() const { return UpperBound; };
+
+    double_t range() const { return UpperBound - LowerBound; };
+
+protected:
+    double_t LowerBound;
+    double_t Nominal;
+    double_t UpperBound;
+};
+
+#endif //OERSTED_COORDINATE_H

src/Optimization/src/CoordinateSpace.cpp

+#include "CoordinateSpace.h"
+
+void CoordinateSpace::limit_velocity(Particle &p) {
+    for (auto &key_coord : Coordinates) {
+        std::string const &key = key_coord.first;
+        Coordinate const &coord = key_coord.second;
+
+        double_t &pv = p.velocity(key);
+        double_t v = SpeedLimit * coord.range();
+        if (pv > v) {
+            pv = v;
+        } else if (pv < -v) {
+            pv = -v;
+        }
+    }
+}
+
+void CoordinateSpace::limit_position(Particle &p) {
+    for (auto &key_coord : Coordinates) {
+        std::string const &key = key_coord.first;
+        Coordinate const &coord = key_coord.second;
+
+        double_t &px = p.position(key);
+        double_t &pv = p.velocity(key);
+
+        double_t lb = coord.lower_bound();
+        double_t ub = coord.upper_bound();
+        if (px < lb) {
+            pv = (lb - px);
+            px = lb;
+        } else if (px > ub) {
+            pv = (ub - px);
+            px = ub;
+        }
+    }
+}
+
+void CoordinateSpace::perturb(Particle &p, std::mt19937 &rng){
+    std::uniform_real_distribution<> dist(-Perturbation,+Perturbation);
+    for (auto &key_coord : Coordinates) {
+        std::string const &key = key_coord.first;
+        Coordinate const &coord = key_coord.second;
+
+        p.position(key) += dist(rng) * coord.range();
+        p.velocity(key) += dist(rng) * coord.range();
+    }
+}
\ No newline at end of file

src/Optimization/src/CoordinateSpace.h

+#ifndef OERSTED_COORDINATESPAC_H
+#define OERSTED_COORDINATESPAC_H
+
+#include "Particle.h"
+#include "Coordinate.h"
+
+class CoordinateSpace {
+public:
+    Coordinate &coordinate(std::string key) { return Coordinates[key]; };
+
+    Particle new_particle(std::mt19937 &rng);
+
+    size_t size() const { return Coordinates.size(); };
+
+    void limit_velocity(Particle &p);
+
+    void limit_position(Particle &p);
+
+    void perturb(Particle &p, std::mt19937 &rng);
+
+    std::map<std::string,Coordinate> Coordinates;
+
+    double_t SpeedLimit{0.2};
+
+    double_t Perturbation{0.2e-3};
+};
+
+
+#endif //OERSTED_COORDINATESPAC_H

src/Optimization/src/CoordinateState.cpp

+#include "CoordinateState.h"

src/Optimization/src/CoordinateState.h

+#ifndef OERSTED_COORDINATESTATE_H
+#define OERSTED_COORDINATESTATE_H
+
+#include <cmath>
+#include <string>
+
+#include "ContiguousDataMap.h"
+
+class CoordinateState {
+public:
+    CoordinateState() {};
+
+    CoordinateState(double_t p, double_t v) : Position{p}, Velocity{v}, PersonalBest{p}, GlobalBest{p} {};
+
+    void update_velocity(double_t omega, double_t phip, double_t phig) { Velocity = omega * (Velocity + phip * (PersonalBest - Position) + phig * (GlobalBest - Position)); }
+
+    void update_position() { Position += Velocity; };
+
+    bool is_converged(double_t tol) { return std::abs(PersonalBest - GlobalBest) < tol; };
+
+    double_t Position;
+
+    double_t Velocity;
+
+    double_t PersonalBest;
+
+    double_t GlobalBest;
+};
+
+
+using ParticleState = ContiguousDataMap<std::string,CoordinateState>;
+
+#endif //OERSTED_COORDINATESTATE_H

src/Optimization/src/Particle.cpp

 #include "Particle.h"
+#include "CoordinateSpace.h"
 
 void Particle::update_personal_best(ObjectiveMap go) {
     update_personal_best(LocalObjective(go));
@@ -90,230 +91,4 @@ Particle CoordinateSpace::new_particle(std::mt19937 &rng) {
     }
 
     return particle;
-}
-
-void CoordinateSpace::limit_velocity(Particle &p) {
-    for (auto &key_coord : Coordinates) {
-        std::string const &key = key_coord.first;
-        Coordinate const &coord = key_coord.second;
-
-        double_t &pv = p.velocity(key);
-        double_t v = SpeedLimit * coord.range();
-        if (pv > v) {
-            pv = v;
-        } else if (pv < -v) {
-            pv = -v;
-        }
-    }
-}
-
-void CoordinateSpace::limit_position(Particle &p) {
-    for (auto &key_coord : Coordinates) {
-        std::string const &key = key_coord.first;
-        Coordinate const &coord = key_coord.second;
-
-        double_t &px = p.position(key);
-        double_t &pv = p.velocity(key);
-
-        double_t lb = coord.lower_bound();
-        double_t ub = coord.upper_bound();
-        if (px < lb) {
-            pv = (lb - px);
-            px = lb;
-        } else if (px > ub) {
-            pv = (ub - px);
-            px = ub;
-        }
-    }
-}
-
-void CoordinateSpace::perturb(Particle &p, std::mt19937 &rng){
-    std::uniform_real_distribution<> dist(-Perturbation,+Perturbation);
-    for (auto &key_coord : Coordinates) {
-        std::string const &key = key_coord.first;
-        Coordinate const &coord = key_coord.second;
-
-        p.position(key) += dist(rng) * coord.range();
-        p.velocity(key) += dist(rng) * coord.range();
-    }
-}
-
-void Swarm::run() {
-    std::uniform_real_distribution<> dist(0.0,1.0);
-
-    bool all_converged{false};
-    for (size_t iter = 0; iter != MaximumIterations && !all_converged; ++iter) {
-        for (size_t i = 0; i != Particles.size(); ++i) {
-            auto go = GlobalObjective(Particles[i]);
-
-            Particles[i].update_personal_best(go);
-            for (size_t j = 0; j != Particles.size(); ++j) {
-                Particles[j].update_global_best(go, Particles[i]);
-            }
-        }
-
-        all_converged = true;
-        for (auto &p : Particles) {
-            all_converged &= p.is_converged(Space,CoordinateTolerance,ObjectiveTolerance);
-        }
-
-        if (!all_converged) {
-            for (auto &p : Particles) {
-                if (!p.has_improved()) {
-                    Space.perturb(p,RNG);
-                }
-            }
-
-            Omega = 2.0 / (1 + M_SQRT2) * dist(RNG);
-
-            for (auto &p : Particles) {
-                p.update_velocity(RNG, Omega, PhiP, PhiG);
-
-                Space.limit_velocity(p);
-
-                p.update_position();
-
-                Space.limit_position(p);
-            }
-        }
-    }
-}
-
-void Swarm::mpirun() {
-    int TAG_PARTICLE_INDEX{1};
-    int TAG_PARTICLE_DATA{2};
-    int TAG_OBJECTIVE_DATA{3};
-    int TAG_EXIT{4};
-    int MPI_RANK;
-    int MPI_SIZE;
-    int FLAG{0};
-
-    MPI_Datatype MPI_CS;
-    MPI_Type_contiguous(sizeof(CoordinateState) / sizeof(MPI_DOUBLE), MPI_DOUBLE, &MPI_CS);
-    MPI_Type_commit(&MPI_CS);
-
-    MPI_Datatype MPI_STATE;
-    MPI_Type_vector((int)Particles[0].state().size(),1,1,MPI_CS,&MPI_STATE);
-    MPI_Type_commit(&MPI_STATE);
-
-
-    MPI_Comm_size(MPI_COMM_WORLD, &MPI_SIZE);
-    MPI_Comm_rank(MPI_COMM_WORLD, &MPI_RANK);
-
-    MPI_Status status;
-    MPI_Request request;
-
-    std::vector<int> PARTICLE_NUMBERS(Particles.size());
-    std::iota(PARTICLE_NUMBERS.begin(),PARTICLE_NUMBERS.end(),0);
-
-    std::deque<int> queue(Particles.size());
-    std::iota(queue.begin(),queue.end(),0);
-
-    MPI_Barrier(MPI_COMM_WORLD);
-
-    if (MPI_RANK == 0) {
-        std::uniform_real_distribution<> dist(0.0,1.0);
-        std::vector<bool> mpi_ready(MPI_SIZE,true);
-        double_t iter{0.0};
-        bool all_ready{false};
-        bool all_converged{false};
-        while ((iter < MaximumIterations && !all_converged) || !all_ready) {
-            // Send particle data to mpi processes
-            for (int i = 1; i != MPI_SIZE && iter < MaximumIterations && !all_converged; ++i) {
-                if (mpi_ready[i] && !queue.empty()) {
-                    int p = queue.front();
-                    queue.pop_front();
-
-                    MPI_Send(PARTICLE_NUMBERS.data() + p, 1, MPI_INT, i, TAG_PARTICLE_INDEX, MPI_COMM_WORLD);
-
-                    MPI_Send(Particles[p].state().vector().data(), 1, MPI_STATE, i, TAG_PARTICLE_DATA, MPI_COMM_WORLD);
-
-                    mpi_ready[i] = false;
-                }
-            }
-
-            // Collect particle data from mpi processes
-            for (int i = 1; i != MPI_SIZE; ++i) {
-                MPI_Iprobe(i, TAG_PARTICLE_INDEX, MPI_COMM_WORLD, &FLAG, &status);
-                if (FLAG) {
-                    int q;
-                    std::vector<double_t> objective(Particles.size(), 0.0);
-
-                    MPI_Recv(&q, 1, MPI_INT, i, TAG_PARTICLE_INDEX, MPI_COMM_WORLD, &status);
-                    MPI_Recv(objective.data(), (int) objective.size(), MPI_DOUBLE, i, TAG_OBJECTIVE_DATA,
-                             MPI_COMM_WORLD, &status);
-
-                    // Update personal best
-                    for (size_t j = 0; j != Particles.size(); ++j) {
-                        Particles[j].update_global_best(objective[j], Particles[q]);
-                    }
-
-                    Particles[q].update_personal_best(objective[q]);
-                    if (!Particles[q].has_improved()) {
-                        Space.perturb(Particles[q],RNG);
-                    }
-
-                    // Update particle q velocity and position
-                    Omega = 2.0 / (1 + M_SQRT2) * dist(RNG);
-                    Particles[q].update_velocity(RNG, Omega, PhiP, PhiG);
-                    Space.limit_velocity(Particles[q]);
-                    Particles[q].update_position();
-                    Space.limit_position(Particles[q]);
-
-                    mpi_ready[i] = true;
-                    queue.push_back(q);
-
-                    iter = iter + 1.0 / Particles.size();
-                }
-            }
-
-            all_converged = true;
-            for (auto &p : Particles) {
-                all_converged &= p.is_converged(Space,CoordinateTolerance,ObjectiveTolerance);
-            }
-
-            all_ready = true;
-            for (bool b : mpi_ready) {
-                all_ready &= b;
-            }
-        }
-
-        // Exit loop
-        for (int i = 1; i != MPI_SIZE; ++i) {
-            MPI_Isend(NULL, 0, MPI_INT, i, TAG_EXIT, MPI_COMM_WORLD, &request);
-        }
-    }
-
-    while (MPI_RANK != 0) {
-        //int data;
-
-        MPI_Iprobe(0,TAG_EXIT,MPI_COMM_WORLD,&FLAG,&status);
-        if (FLAG) {
-            MPI_Irecv(NULL,0,MPI_INT,0,TAG_EXIT,MPI_COMM_WORLD,&request);
-            break;
-        }
-
-        MPI_Iprobe(0,TAG_PARTICLE_INDEX,MPI_COMM_WORLD,&FLAG,&status);
-        if (FLAG) {
-            int p;
-
-            MPI_Recv(&p,1,MPI_INT,0,TAG_PARTICLE_INDEX,MPI_COMM_WORLD,&status);
-            MPI_Recv(Particles[p].state().vector().data(),1,MPI_STATE,0,TAG_PARTICLE_DATA,MPI_COMM_WORLD,&status);
-
-            Particles[p].state(0).PersonalBest = MPI_RANK;
-
-            ObjectiveMap go = GlobalObjective(Particles[p]);
-
-            std::vector<double_t> objective(Particles.size(),0.0);
-            for (size_t i = 0; i != Particles.size(); ++i) {
-                objective[i] = Particles[i].local_objective(go);
-            }
-
-            MPI_Send(&p,1,MPI_INT,0,TAG_PARTICLE_INDEX,MPI_COMM_WORLD);
-            MPI_Send(objective.data(),(int)objective.size(),MPI_DOUBLE,0,TAG_OBJECTIVE_DATA,MPI_COMM_WORLD);
-        }
-    }
-
-    MPI_Type_free(&MPI_CS);
-    MPI_Type_free(&MPI_STATE);
 }
\ No newline at end of file

src/Optimization/src/Particle.h

@@ -12,84 +12,10 @@
 #include <mpi.h>
 
 
-class Particle;
-class CoordinateSpace;
-class Coordinate;
-class CoordinateState;
-
-template <class Key, class T>
-class ContiguousDataMap {
-public:
-    // TODO: Copy Constructor
-
-    T &operator[] (Key k) {
-        auto it = Map.find(k);
-        if (it != Map.end()) {
-            return Vector[it->second];
-        } else {
-            Map[k] = Vector.size();
-            Vector.push_back(T());
-            return Vector.back();
-        }
-    };
-
-    T &operator[] (size_t k) { return Vector[k]; };
-
-    std::map<Key,size_t> &map() { return Map; };
-
-    std::vector<T> &vector() { return Vector; };
-
-    size_t size() const { return Vector.size(); };
-
-protected:
-    std::map<Key,size_t> Map;
-
-    std::vector<T> Vector;
-};
-
-using ObjectiveMap = ContiguousDataMap<std::string,double_t>;
-using ParticleState = ContiguousDataMap<std::string,CoordinateState>;
-
-class CoordinateState {
-public:
-    CoordinateState() {};
-
-    CoordinateState(double_t p, double_t v) : Position{p}, Velocity{v}, PersonalBest{p}, GlobalBest{p} {};
-
-    void update_velocity(double_t omega, double_t phip, double_t phig) { Velocity = omega * (Velocity + phip * (PersonalBest - Position) + phig * (GlobalBest - Position)); }
-
-    void update_position() { Position += Velocity; };
-
-    bool is_converged(double_t tol) { return std::abs(PersonalBest - GlobalBest) < tol; };
-
-    double_t Position;
-
-    double_t Velocity;
-
-    double_t PersonalBest;
-
-    double_t GlobalBest;
-};
-
-class Coordinate {
-public:
-    Coordinate() {};
-
-    Coordinate(double_t lb, double_t nom, double_t ub) : LowerBound{lb}, Nominal{nom}, UpperBound{ub} {};
-
-    double_t lower_bound() const { return LowerBound; };
-
-    double_t nominal() const { return Nominal; };
-
-    double_t upper_bound() const { return UpperBound; };
+#include "CoordinateState.h"
+#include "ContiguousDataMap.h"
 
-    double_t range() const { return UpperBound - LowerBound; };
-
-protected:
-    double_t LowerBound;
-    double_t Nominal;
-    double_t UpperBound;
-};
+class CoordinateSpace;
 
 class Particle {
 public:
@@ -146,61 +72,4 @@ protected:
 
 };
 
-
-class CoordinateSpace {
-public:
-    Coordinate &coordinate(std::string key) { return Coordinates[key]; };
-
-    Particle new_particle(std::mt19937 &rng);
-
-    size_t size() const { return Coordinates.size(); };
-
-    void limit_velocity(Particle &p);
-
-    void limit_position(Particle &p);
-
-    void perturb(Particle &p, std::mt19937 &rng);
-
-    std::map<std::string,Coordinate> Coordinates;
-
-    double_t SpeedLimit{0.2};
-
-    double_t Perturbation{0.2e-3};
-};
-
-class Swarm {
-public:
-    Swarm(CoordinateSpace ss, std::function<ObjectiveMap(Particle)> go, size_t Np, double_t otol, size_t maxit) : Space{ss}, GlobalObjective{go}, ObjectiveTolerance{otol}, MaximumIterations{maxit} {
-        std::random_device rd;
-        RNG = std::mt19937{rd()};
-
-        for (size_t i = 0; i != Np; ++i) {
-            Particles.push_back(Space.new_particle(RNG));
-        }
-    }
-
-    size_t size() const { return Particles.size(); };
-
-    void run();
-
-    void mpirun();
-
-    std::vector<Particle> Particles;
-
-    double_t Omega{1.0 / (1.0 + M_SQRT2)};
-    double_t PhiP{4.0};
-    double_t PhiG{4.0};
-
-    double_t CoordinateTolerance{1e-2};
-    double_t ObjectiveTolerance;
-    size_t MaximumIterations;
-
-protected:
-    CoordinateSpace Space;
-
-    std::function<ObjectiveMap(Particle)> GlobalObjective;
-
-    std::mt19937 RNG;
-};
-
 #endif //OERSTED_PARTICLE_H

src/Optimization/src/Swarm.cpp

+#include "Swarm.h"
+
+void Swarm::run() {
+    std::uniform_real_distribution<> dist(0.0,1.0);
+
+    bool all_converged{false};
+    for (size_t iter = 0; iter != MaximumIterations && !all_converged; ++iter) {
+        for (size_t i = 0; i != Particles.size(); ++i) {
+            auto go = GlobalObjective(Particles[i]);
+
+            Particles[i].update_personal_best(go);
+            for (size_t j = 0; j != Particles.size(); ++j) {
+                Particles[j].update_global_best(go, Particles[i]);
+            }
+        }
+
+        all_converged = true;
+        for (auto &p : Particles) {
+            all_converged &= p.is_converged(Space,CoordinateTolerance,ObjectiveTolerance);
+        }
+
+        if (!all_converged) {
+            for (auto &p : Particles) {
+                if (!p.has_improved()) {
+                    Space.perturb(p,RNG);
+                }
+            }