Skip to content

Commit

Permalink
Unshadow members; Explicit casts; Tidy tests
Browse files Browse the repository at this point in the history
  • Loading branch information
pmslavin committed Dec 31, 2024
1 parent d378b23 commit 8b7dc43
Show file tree
Hide file tree
Showing 7 changed files with 130 additions and 100 deletions.
28 changes: 14 additions & 14 deletions include/pagmo/algorithms/nsga3.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -38,23 +38,23 @@ class PAGMO_DLL_PUBLIC nsga3{
const log_type &get_log() const { return m_log; }
void set_verbosity(unsigned level) { m_verbosity = level; }
unsigned get_verbosity() const { return m_verbosity; }
void set_seed(unsigned seed) { reng.seed(seed); seed = seed; }
unsigned get_seed() const { return seed; }
bool has_memory() const {return use_memory; }
void set_seed(unsigned seed) { m_reng.seed(seed); m_seed = seed; }
unsigned get_seed() const { return m_seed; }
bool has_memory() const {return m_use_memory; }
private:
unsigned ngen;
double cr; // crossover
double eta_c; // eta crossover
double mut; // mutation
double eta_mut; // eta mutation
size_t divisions; // Reference Point hyperplane subdivisions
unsigned seed; // Seed for PRNG initialisation
bool use_memory; // Preserve extremes, ideal, nadir across generations
mutable NSGA3Memory memory{};
mutable detail::random_engine_type reng; // Defaults to std::mt19937
unsigned m_gen;
double m_cr; // crossover
double m_eta_c; // eta crossover
double m_mut; // mutation
double m_eta_mut; // eta mutation
size_t m_divisions; // Reference Point hyperplane subdivisions
unsigned m_seed; // Seed for PRNG initialisation
bool m_use_memory; // Preserve extremes, ideal, nadir across generations
mutable NSGA3Memory m_memory{};
mutable detail::random_engine_type m_reng; // Defaults to std::mt19937
mutable log_type m_log;
unsigned m_verbosity;
mutable std::vector<ReferencePoint> refpoints;
mutable std::vector<ReferencePoint> m_refpoints;
// Serialisation support
friend class boost::serialization::access;
template <typename Archive>
Expand Down
1 change: 0 additions & 1 deletion include/pagmo/utils/multi_objective.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -227,7 +227,6 @@ PAGMO_DLL_PUBLIC double perpendicular_distance(const std::vector<double> &, cons

/* Choose single random element from vector container */
template <class T>
PAGMO_DLL_PUBLIC
T choose_random_element(const std::vector<T> &container){
std::vector<T> choice;
std::sample(container.begin(), container.end(), std::back_inserter(choice), 1, std::mt19937{random_device::next()});
Expand Down
6 changes: 3 additions & 3 deletions include/pagmo/utils/reference_point.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -21,9 +21,9 @@ namespace pagmo{
class PAGMO_DLL_PUBLIC ReferencePoint{
public:
ReferencePoint(size_t nobj);
~ReferencePoint();
size_t dim() const;
double& operator[](int);
//double& operator[](int);
double& operator[](size_t);
friend PAGMO_DLL_PUBLIC std::ostream& operator<<(std::ostream& ostr, const ReferencePoint& rp);
void increment_members(){ ++nmembers; }
void decrement_members(){ --nmembers; }
Expand Down Expand Up @@ -56,7 +56,7 @@ void associate_with_reference_points(

size_t identify_niche_point(std::vector<ReferencePoint> &);

size_t n_choose_k(unsigned, unsigned);
size_t n_choose_k(size_t, size_t);

} // namespace pagmo

Expand Down
78 changes: 38 additions & 40 deletions src/algorithms/nsga3.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -28,8 +28,8 @@ namespace pagmo{

nsga3::nsga3(unsigned gen, double cr, double eta_c, double mut, double eta_mut,
size_t divisions, unsigned seed, bool use_memory)
: ngen(gen), cr(cr), eta_c(eta_c), mut(mut), eta_mut(eta_mut),
divisions(divisions), seed(seed), use_memory(use_memory), reng(seed){
: m_gen(gen), m_cr(cr), m_eta_c(eta_c), m_mut(mut), m_eta_mut(eta_mut),
m_divisions(divisions), m_seed(seed), m_use_memory(use_memory), m_reng(seed){
// Validate ctor args
if(cr < 0.0 || cr > 1.0){
pagmo_throw(std::invalid_argument, "The crossover probability must be in the range [0, 1], while a value of "
Expand Down Expand Up @@ -63,11 +63,11 @@ nsga3::nsga3(unsigned gen, double cr, double eta_c, double mut, double eta_mut,

std::vector<ReferencePoint> nsga3::generate_uniform_reference_points(size_t nobjs, size_t divisions) const{
ReferencePoint rp(nobjs);
if(!refpoints.empty()){
refpoints.clear();
if(!m_refpoints.empty()){
m_refpoints.clear();
}
refpoints = generate_reference_point_level(rp, divisions, 0, divisions);
return refpoints;
m_refpoints = generate_reference_point_level(rp, divisions, 0, divisions);
return m_refpoints;
}


Expand All @@ -77,14 +77,14 @@ std::vector<std::vector<double>> nsga3::translate_objectives(population pop) con
auto objs = pop.get_f();
std::vector<double> p_ideal;
if(has_memory()){
decltype(objs) combined{memory.v_ideal};
if(memory.v_ideal.size() != 0){ // i.e. not first gen
decltype(objs) combined{m_memory.v_ideal};
if(m_memory.v_ideal.size() != 0){ // i.e. not first gen
combined.insert(combined.end(), objs.begin(), objs.end());
p_ideal = ideal(combined);
}else{
p_ideal = ideal(objs);
}
memory.v_ideal = p_ideal;
m_memory.v_ideal = p_ideal;
}else{
p_ideal = ideal(objs);
}
Expand Down Expand Up @@ -113,16 +113,16 @@ std::vector<std::vector<double>> nsga3::find_extreme_points(population pop,
std::vector<double> min_obj{};

if(has_memory()){
if(memory.v_extreme.size() == 0){
if(m_memory.v_extreme.size() == 0){
for(size_t idx=0; idx<nobj ; idx++){
memory.v_extreme.push_back(std::vector<double>(nobj, {}));
m_memory.v_extreme.push_back(std::vector<double>(nobj, {}));
}
}else{
for(size_t p=0; p<memory.v_extreme.size(); p++){
double asf = achievement(memory.v_extreme[p], weights);
for(size_t p=0; p<m_memory.v_extreme.size(); p++){
double asf = achievement(m_memory.v_extreme[p], weights);
if(asf < min_asf){
min_asf = asf;
min_obj = memory.v_extreme[p];
min_obj = m_memory.v_extreme[p];
}
}
}
Expand All @@ -139,7 +139,7 @@ std::vector<std::vector<double>> nsga3::find_extreme_points(population pop,
}
points.push_back(std::vector<double>(min_obj));
if(has_memory()){
memory.v_extreme[i] = std::vector<double>(min_obj);
m_memory.v_extreme[i] = std::vector<double>(min_obj);
}
}

Expand Down Expand Up @@ -195,14 +195,14 @@ std::vector<double> nsga3::find_intercepts(population pop, std::vector<std::vect
auto objs = pop.get_f();
std::vector<double> v_nadir;
if(has_memory()){
decltype(objs) combined{memory.v_nadir};
if(memory.v_nadir.size() != 0){
decltype(objs) combined{m_memory.v_nadir};
if(m_memory.v_nadir.size() != 0){
combined.insert(combined.end(), objs.begin(), objs.end());
v_nadir = nadir(combined);
}else{
v_nadir = nadir(objs);
}
memory.v_nadir = v_nadir;
m_memory.v_nadir = v_nadir;
}else{
v_nadir = nadir(objs);
}
Expand Down Expand Up @@ -241,8 +241,6 @@ population nsga3::evolve(population pop) const{
auto dim_i = prob.get_nix();
auto NP = pop.size();

// Initialize the population

/* Verify problem characteristics:
* - Has multiple objectives
* - Is not stochastic
Expand Down Expand Up @@ -271,11 +269,11 @@ population nsga3::evolve(population pop) const{
"NSGA-III requires a population greater than 5 and which is divisible by 4."
"Detected input population size is: " + std::to_string(NP));
}
size_t num_rps = n_choose_k(prob.get_nf() + divisions - 1, divisions);
size_t num_rps = n_choose_k(prob.get_nf() + m_divisions - 1, m_divisions);
if(NP <= num_rps){
pagmo_throw(std::invalid_argument,
"Population size must exceed number of reference points. NP = "
+ std::to_string(NP) + " while " + std::to_string(divisions) + " divisions for "
+ std::to_string(NP) + " while " + std::to_string(m_divisions) + " divisions for "
"reference points gives a total of " + std::to_string(num_rps) + " points.");
}

Expand All @@ -289,13 +287,13 @@ population nsga3::evolve(population pop) const{
std::iota(shuffle1.begin(), shuffle1.end(), vector_double::size_type(0));
std::iota(shuffle2.begin(), shuffle2.end(), vector_double::size_type(0));

for(decltype(ngen)gen = 1u; gen <= ngen; gen++){
for(decltype(m_gen)gen = 1u; gen <= m_gen; gen++){
// Copy existing population
population popnew(pop);

// Permute population indices
std::shuffle(shuffle1.begin(), shuffle1.end(), reng);
std::shuffle(shuffle2.begin(), shuffle2.end(), reng);
std::shuffle(shuffle1.begin(), shuffle1.end(), m_reng);
std::shuffle(shuffle2.begin(), shuffle2.end(), m_reng);

/* 1. Generate offspring population Q_t
* 2. R = P_t U Q_t
Expand All @@ -304,7 +302,7 @@ population nsga3::evolve(population pop) const{

if(m_verbosity > 0u){
std::vector<double> p_ideal = ideal(pop.get_f());
if (gen % m_verbosity == 1u || m_verbosity == 1u) {
if (m_gen % m_verbosity == 1u || m_verbosity == 1u) {
// We compute the ideal point
// Every 50 lines print the column names
if (count % 50u == 1u) {
Expand All @@ -318,7 +316,7 @@ population nsga3::evolve(population pop) const{
}
print('\n');
}
print(std::setw(7), gen, std::setw(15), prob.get_fevals() - fevals0);
print(std::setw(7), m_gen, std::setw(15), prob.get_fevals() - fevals0);
for (decltype(p_ideal.size()) i = 0u; i < p_ideal.size(); ++i) {
if (i >= 5u) {
break;
Expand All @@ -328,18 +326,18 @@ population nsga3::evolve(population pop) const{
print('\n');
++count;
}
m_log.emplace_back(gen, prob.get_fevals() - fevals0, p_ideal);
m_log.emplace_back(m_gen, prob.get_fevals() - fevals0, p_ideal);
}

// Offspring generation
for (decltype(NP) i = 0; i < NP; i += 4) {
// We create two offsprings using the shuffled list 1
decltype(shuffle1) parents1;
std::sample(shuffle1.begin(), shuffle1.end(), std::back_inserter(parents1), 2, std::mt19937{reng()});
std::sample(shuffle1.begin(), shuffle1.end(), std::back_inserter(parents1), 2, std::mt19937{m_reng()});
children = detail::sbx_crossover_impl(pop.get_x()[parents1[0]], pop.get_x()[parents1[1]], bounds, dim_i,
cr, eta_c, reng);
detail::polynomial_mutation_impl(children.first, bounds, dim_i, mut, eta_mut, reng);
detail::polynomial_mutation_impl(children.second, bounds, dim_i, mut, eta_mut, reng);
m_cr, m_eta_c, m_reng);
detail::polynomial_mutation_impl(children.first, bounds, dim_i, m_mut, m_eta_mut, m_reng);
detail::polynomial_mutation_impl(children.second, bounds, dim_i, m_mut, m_eta_mut, m_reng);
// Evaluation via prob ensures feval counter is correctly updated
auto f1 = prob.fitness(children.first);
auto f2 = prob.fitness(children.second);
Expand All @@ -348,11 +346,11 @@ population nsga3::evolve(population pop) const{

// Repeat with the shuffled list 2
decltype(shuffle2) parents2;
std::sample(shuffle2.begin(), shuffle2.end(), std::back_inserter(parents2), 2, std::mt19937{reng()});
std::sample(shuffle2.begin(), shuffle2.end(), std::back_inserter(parents2), 2, std::mt19937{m_reng()});
children = detail::sbx_crossover_impl(pop.get_x()[parents2[0]], pop.get_x()[parents2[1]], bounds, dim_i,
cr, eta_c, reng);
detail::polynomial_mutation_impl(children.first, bounds, dim_i, mut, eta_mut, reng);
detail::polynomial_mutation_impl(children.second, bounds, dim_i, mut, eta_mut, reng);
m_cr, m_eta_c, m_reng);
detail::polynomial_mutation_impl(children.first, bounds, dim_i, m_mut, m_eta_mut, m_reng);
detail::polynomial_mutation_impl(children.second, bounds, dim_i, m_mut, m_eta_mut, m_reng);
f1 = prob.fitness(children.first);
f2 = prob.fitness(children.second);
popnew.push_back(children.first, f1);
Expand Down Expand Up @@ -389,7 +387,7 @@ std::vector<size_t> nsga3::selection(population &R, size_t N_pop) const{
while(next_size < N_pop){
next_size += fronts[last_front++].size();
}
fronts.erase(fronts.begin() + last_front, fronts.end());
fronts.erase(fronts.begin() + static_cast<std::vector<vector_double>::difference_type>(last_front), fronts.end());

// Accept all members of first l-1 fronts
for(size_t f=0; f<fronts.size()-1; f++){
Expand All @@ -406,7 +404,7 @@ std::vector<size_t> nsga3::selection(population &R, size_t N_pop) const{
auto ext_points = find_extreme_points(R, fronts, translated_objectives);
auto intercepts = find_intercepts(R, ext_points);
auto norm_objs = normalize_objectives(translated_objectives, intercepts);
std::vector<ReferencePoint> rps = generate_uniform_reference_points(nobj, divisions);
std::vector<ReferencePoint> rps = generate_uniform_reference_points(nobj, m_divisions);
associate_with_reference_points(rps, norm_objs, fronts);

// Apply RP selection to final front until N_pop reached
Expand All @@ -418,7 +416,7 @@ std::vector<size_t> nsga3::selection(population &R, size_t N_pop) const{
rps[min_rp_idx].remove_candidate(selected_idx.value());
next.push_back(selected_idx.value());
}else{
rps.erase(rps.begin() + min_rp_idx);
rps.erase(rps.begin() + static_cast<std::vector<vector_double>::difference_type>(min_rp_idx));
}
}

Expand All @@ -428,7 +426,7 @@ std::vector<size_t> nsga3::selection(population &R, size_t N_pop) const{
// Object serialization
template <typename Archive>
void nsga3::serialize(Archive &ar, unsigned int) {
detail::archive(ar, ngen, cr, eta_c, mut, eta_mut, seed, m_verbosity, m_log);
detail::archive(ar, m_gen, m_cr, m_eta_c, m_mut, m_eta_mut, m_seed, m_verbosity, m_log);
}

} // namespace pagmo
Expand Down
40 changes: 36 additions & 4 deletions src/utils/multi_objective.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -638,7 +638,16 @@ vector_double decompose_objectives(const vector_double &f, const vector_double &
return {fd};
}

/* Gaussian Elimination */
/// Gaussian Elimination
/**
* Solves the linear system Ax = b by performing Gaussian elimination
* on an augmented matrix A|b.
*
* @param
* @param
* @return
* @throws
*/
vector_double gaussian_elimination(std::vector<std::vector<double>> A, const vector_double &b){
// N.B. Check dimensions and pivots
if(A[0][0] == 0.0){
Expand All @@ -661,18 +670,32 @@ vector_double gaussian_elimination(std::vector<std::vector<double>> A, const vec
}
vector_double x(N);
// Back substitution
for(int i=N-1; i>=0; i--){
size_t i = N-1;
while(true){
for(size_t var=i+1; var<N; var++){
A[i][N] -= A[i][var]*x[var];
}
x[i] = A[i][N]/A[i][i];
if(i == 0){
break;
}
i--;
}

return x;
}


/* Achievement Scalarization Function */
/// Achievement Scalarization Function
/**
*
*
*
* @param
* @param
* @return
* @throws
*/
double achievement(const vector_double &objs, const vector_double &weights){
const double default_weight = 1e-5;
double max_ratio = -std::numeric_limits<double>::max();
Expand All @@ -686,7 +709,16 @@ double achievement(const vector_double &objs, const vector_double &weights){
return max_ratio;
}

/* Distance from objective point to perpendicular intersection with reference point vector */
/// Distance from objective point to perpendicular intersection with reference point vector */
/**
*
*
*
* @param
* @param
* @return
* @throws
*/
double perpendicular_distance(const std::vector<double> &ref_point, const std::vector<double> &obj_point){
double num = 0.0, denom = 0.0, sq_dist = 0.0;
for(size_t i=0; i<ref_point.size(); i++){
Expand Down
Loading

0 comments on commit 8b7dc43

Please sign in to comment.