Эффективные алгоритмы инкрементальной и асинхронной недоминирующей сортировки для многокритериальных эволюционных алгоритмов тема диссертации и автореферата по ВАК РФ 05.13.11, кандидат наук Якупов Илья Юрьевич

• Якупов Илья Юрьевич
• кандидат наук
• 2019

Оглавление диссертации кандидат наук Якупов Илья Юрьевич

Реферат

Synopsis

Introduction

Chapter 1. Analytical Overview and Motivation

1.1 Mathematical Optimization

1.1.1 Problem Overview

1.1.2 Approaches for Multiobjective Optimization Problems

1.2 Notation and Basic Definitions

1.3 Generational Multiobjective Optimization Algorithms

1.3.1 Aggregation-Based Algorithms

1.3.2 Criterion-Based Algorithms

1.3.3 Dominance-Based Algorithms

1.3.4 Hybrid Algorithms

1.3.5 The NSGA-II Algorithm

1.4 Steady-state evolutionary multiobjective optimization algorithms

1.4.1 Incremental Non-Dominated Sorting

1.4.2 Efficient Non-domination Level Update

1.5 Chapter 1 Conclusion

Chapter 2. Two-Dimensional Incremental Non-Dominated Sorting

2.1 Used Data Structures

2.2 Algorithm Description

2.2.1 Data Structure

2.2.2 Lookup

2.2.3 Insertion

2.2.4 Deletion of the Worst Solution

2.2.5 Querying the Л-th Solution

2.3 Experiments

2.3.1 Synthetic Tests

2.3.2 Tests as a Part of NSGA-II

2.3.3 Scalability

2.4 Chapter 2 Conclusion

Chapter 3. Incremental Non-Dominated Sorting for Arbitrary Dimensions

3.1 Preliminaries

3.1.1 Definitions

3.1.2 The Divide-and-Conquer Approach

3.2 Adapting the Divide-and-Conquer Algorithm to Incremental Sorting

3.2.1 Finding the Rank of the New Point

3.2.2 Reducing the Set of Points to Work With

3.2.3 Taking Advantage of Known Ranks

3.2.4 Is-Anything-Changed Heuristic

3.2.5 Putting Everything Together

3.3 Experiments

3.4 Chapter 3 Conclusion

Chapter 4. Algorithms for Asynchronous Non-Dominated Sorting

4.1 Concurrency Primitives

4.2 Introducing Concurrency

4.2.1 The Compare-And-Set Approach

4.2.2 A Time-Stamping Modification

4.2.3 The Approach with Finer-Grained Locks

4.2.4 Recomputation of the Crowding Distance

4.3 Experiments

4.3.1 Asynchronous NSGA-II: Details and Parameters

4.3.2 Experiments with Fixed Number of Threads

4.3.3 Results and Discussion: Running Times

4.3.4 Results and Discussion: Hypervolumes

4.3.5 Exploring Scalability When Adding Worker Threads

4.3.6 Performance under Heavy CPU Load

4.4 Chapter 4 Conclusion

Conclusion

References

List of Figures

List of Tables

Appendix A. Copies of Author's Publications

Реферат

Введение диссертации (часть автореферата) на тему «Эффективные алгоритмы инкрементальной и асинхронной недоминирующей сортировки для многокритериальных эволюционных алгоритмов»

Общая характеристика работы

Актуальность темы. Задача многокритериальной оптимизации состоит в одновременной оптимизации нескольких, возможно, конфликтующих, целевых функций в некоторой области их определения. А именно, в данной задаче требуется найти вектор значений целевых переменных, удовлетворяющий наложенным ограничениям и оптимизирующий векторную функцию, элементы которой соответствуют целевым функциям. Отметим, что многие задачи оптимизации, встречающиеся на практике, являются многокритериальными, и их решение в каких-либо однокритериальных формулировках не дает полноты знаний о возможных решениях задачи.

Критерии оптимальности выбираются в зависимости от поставленной задачи. Одним из наиболее часто рассматриваемых критериев является оптимальность по Парето. Решение, оптимальное по Парето, обладает следующим свойством: значение каждой целевой функции не может быть улучшено без ухудшения значения по крайней мере одной другой функции. Множество всех решений задачи, оптимальных по Парето, называется множеством Парето, или же множеством неулучшаемых альтернатив. Множество же векторов значений целевых функций, соответствующих множеству Парето, называют фронтом Парето.

Без ограничения общности будем считать, что каждую целевую функцию требуется минимизировать, а значениями каждой целевой функции являются вещественные числа.

Важной концепцией является доминирование по Парето. В предположении минимизации целевых функций оно определяется следующим образом. Решение задачи многокритериальной оптимизации, имеющее вектор значений целевых функций р, строго доминирует по Парето решение с вектором ц, что обозначается как р ■ ц, тогда и только тогда, когда для каждой целевой функции г ее значение для первого из решений рг не превосходит ее же значения для второго решения цг, и существует хотя бы одна целевая функция, для которой это неравенство выполняется строго:

14

1 < г < к,рг < Цг] 1 < г < к,рг < цг.

Таким образом, множество Парето для задачи оптимизации можно также определить как множество решений этой задачи, каждое из которых не доминируется ни одним другим ее решением.

И множество Парето, и фронт Парето практически всех многокритериальных задач оптимизации содержит более одного решения ввиду конфликтов между целевыми функциями. Так как знание какого-либо одного решения, принадлежащего множеству Парето, недостаточно для того, чтобы сделать выводы о всех возможных оптимальных решениях этой задачи, в многокритериальной постановке, как правило, требуется найти множество решений, покрывающих, по возможности равномерно, весь фронт Парето. В связи с этим, в многокритериальной оптимизации широкое применение находят эволюционные алгоритмы, так как они работают одновременно с множеством решений задачи. В терминологии эволюционных алгоритмов, решение задачи называется особью, совокупность одновременно рассматриваемых решений — популяцией, целевая функция — функцией приспособленности1.

Современные многокритериальные эволюционные алгоритмы можно разделить на три класса: алгоритмы, напрямую основанные на отношении доминирования по Парето (примеры: NSGA-II2, SPEA23, PESA-II4); алгоритмы, оптимизирующие так называемые индикаторы — функции, характеризующие одним числом качество популяции в целом (примеры: IBEA5, HypE6); алгоритмы, сводящие многокритериальную задачу к нескольким скалярным задачам (пример: MOEA/D7); также известны гибридные алгоритмы, такие какNSGA-Ш8.

Значительная часть алгоритмов, относящихся к первому из указанных классов (а также гибридные алгоритмы), используют для ранжирования решений процедуру, известную как недоминирующая сортировка9 (англ. non-dominated sorting). Данная процедура использует только векторы целевых функций решений задачи. Пусть дано множество (вообще говоря, мультимножество) решений задачи P. Каждому из них назначается ранг по следующему правилу:

- каждому решению p, такому что yq G P выполняется условие q ^ p, назначается нулевой ранг: r(p) ^ 0;

lLuke, S. Essentials of Metaheuristics. Lulu, 2009. 253 p.

2Deb, K. [et al.]. A Fast and Elitist Multi-Objective Genetic Algorithm: NSGA-II // IEEE Transactions on Evolutionary Computation. 2002. Vol. 6, no. 2. P. 182-197.

3Zitzler, E., Laumanns, M., Thiele, L. SPEA2: Improving the Strength Pareto Evolutionary Algorithm for Multiobjective Optimization//Proceedings of the EUR0GEN'2001 Conference. 2001. P. 95-100.

4Corne, D. W. [etal.]. PESA-II: Region-based Selection in Evolutionary Multiobjective Optimization // Proceedings of Genetic and Evolutionary Computation Conference. 2001. P. 283-290.

5Zitzler, E., Kunzli, S. Indicator-Based Selection in Multiobjective Search // Parallel Problem Solving from Nature - PPSN VIII. 2004. P. 832-842. (Lecture Notes in Computer Science ; 3242).

6Bader, J., Zitzler, E. HypE: An Algorithm for Fast Hypervolume-Based Many-Objective Optimization//Evolutionary Computation. 2011. Vol. 19, no. 1. P. 45-76.

1 Zhang, Q., Li, H. MOEA/D: A Multiobjective Evolutionary Algorithm Based on Decomposition // IEEE Transactions on Evolutionary Computation. 2007. Vol. 11, no. 6. P. 712-731.

8Deb, K., Jain, H. An Evolutionary Many-Objective Optimization Algorithm Using Reference-Point-Based Nondominated Sorting Approach, Part I: Solving Problems With Box Constraints // IEEE Transactions on Evolutionary Computation. 2013. Vol. 18, no. 4. P. 577-601.

9 Имеется также вариант перевода «недоминируемая сортировка». Оба варианта перевода, как, впрочем, и исходное название на английском языке, не обладают высокими литературными качествами, но являются достаточно краткими и удобными в употреблении.

- для всех остальных решений ранг назначается следующим образом: r(p) ^ 1 + maxgePig^p r(q).

Таким образом, чем больше ранг решения, тем длиннее цепочка доминирования от него до некоторого решения задачи, не доминируемого другими решениями из данного множества, и тем менее качественным признается данное решение. Процедура недоминирующей сортировки обладает относительно высокой вычислительной сложностью: большая часть известных для нее алгоритмов имеют сложность ©(n2k), где n — число решений, к — число целевых функций.

Эволюционные алгоритмы имеют репутацию алгоритмов оптимизации, хорошо поддающихся распараллеливанию. Однако алгоритмы, построенные по классической, популяционной схеме, состоящей в синтезе большого числа решений и вычислении их приспособленности, вынуждены ожидать окончания вычисления приспособленности всех новых решений. Если время вычисления приспособленности не является постоянным, данное обстоятельство существенно снижает коэффициент использования ресурсов многопроцессорных и распределенных вычислительных систем. Для преодоления данного недостатка были предложены инкрементальные, или стационарные эволюционные алгоритмы (англ. steady-state evolutionary algorithms), в которых на каждой итерации синтезируется одна (или небольшое число) особей, вычисляется их приспособленность, а затем обновляется состояние алгоритма. Алгоритмы указанного типа могут быть преобразованы в алгоритмы с асинхронным вычислением приспособленности: в каждый момент времени в состоянии вычисления приспособленности может находиться несколько решений, каждое из которых обрабатывается в отдельном потоке. Решения с вычисленной приспособленностью интегрируется в состояние алгоритма сразу же после вычисления, без необходимости ожидания вычисления приспособленности остальных особей. Отметим также, что инкрементальные эволюционные алгоритмы нередко имеют преимущество и в качестве получаемых решений10.

Главным недостатком инкрементальных (и, следовательно, асинхронных) эволюционных алгоритмов является то, что обновление состояния алгоритма, выполняемое при поступлении каждого решения с вычисленной приспособленностью, может быть более длительным в пересчете на одну особь. Вычислительно сложные процедуры обновления состояния алгоритма, такие как недоминирующая сортировка, характерны для многокритериальных алгоритмов. Поэтому тема настоящей диссертации, посвященной алгоритмам и структурам данных для эффективной недоминирующей сортировке в условиях инкрементальных и асинхронных обновлений, является актуальной.

Степень разработанности темы. В свете того, что тема диссертации не затрагивает вопросов о целевых функциях, а следовательно, о соответствии векторов значений целевых переменных («генотипа» особи) и векторов значений

10 Nebro, A. J., Durillo, J. J. On the Effect of Applying a Steady-State Selection Scheme in the Multi-Objective Genetic Algorithm NSGA-II // Nature-Inspired Algorithms for Optimisation. 2009. P. 435-456. (Studies in Computational Intelligence ; 193).

целевых функций (приспособленности особи), далее будем отождествлять решение (особь) и его приспособленность (и рассматривать только приспособленность), что позволяет считать решения точками в к-мерном пространстве значений целевых функций.

Процедура недоминирующей сортировки была предложена для использования в многокритериальных эволюционных алгоритмах в качестве части алгоритма NSGA11, в этой же работе был предложен алгоритм со сложностью <Э(п3к) по времени и О(п) по памяти, где п — число точек. В следующей версии этого алгоритма, NSGA-П12, алгоритм для недоминирующей сортировки был улучшен, сложность его составила 0(п2к) по времени и О(п2) по памяти, что отчасти способствовало популярности алгоритма NSGA-П. Однако указанная сложность остается высокой; так, практически во всех алгоритмах, использующих недоминирующую сортировку, она остается либо самым узким местом алгоритма по асимптотике, либо одним из таких мест.

В связи с этим в ряде последующих работ были предложены альтернативные, более эффективные алгоритмы для недоминирующей сортировки. Большая часть из них имеет сложность <Э(п2к) по времени в худшем случае и О(пк) по памяти13. Однако алгоритм, предложенный в работе Йенсена14 и впоследствии улучшенный в ряде недавних работ15, имеет сложность по времени O(n(logп)к-1) и большую, по сравнению с другими алгоритмами, эффективность при большом числе точек.

Вопрос об эффективности недоминирующей сортировки как части инкрементального многокритериального эволюционного алгоритма был поднят в ра-

nSrinivas, N., Deb, K. Multiobjective Optimization Using Nondominated Sorting in Genetic Algorithms // Evolutionary Computation. 1994. Vol. 2, no. 3. P. 221-248.

12Deb, K. [et al.]. A Fast and Elitist Multi-Objective Genetic Algorithm: NSGA-II // IEEE Transactions on Evolutionary Computation. 2002. Vol. 6, no. 2. P. 182-197.

13Roy, P. C., Islam, M. M., Deb, K. Best Order Sort: A New Algorithm to Non-dominated Sorting for Evolutionary Multi-objective Optimization // Proceedings of Genetic and Evolutionary Computation Conference Companion. 2016. P. 1113-1120; Roy, P. C, Deb, K., Islam, M.M. An Efficient Nondominated Sorting Algorithm for Large Number of Fronts // IEEE Transactions on Cybernetics. 2019. Vol. 49, no. 3. P. 859-869; Gustavsson, P., Syberfeldt, A. A New Algorithm Using the Non-dominated Tree to Improve Non-dominated Sorting//Evolutionary Computation. 2018. Vol. 26, no. 1. P. 89-116; Zhang, X.[et al.]. An Efficient Approach to Nondominated Sorting for Evolutionary Multiobjective Optimization // IEEE Transactions on Evolutionary Computation. 2015. Vol. 19, no. 2. P. 201-213.

14Jensen, M. T. Reducing the Run-time Complexity of Multiobjective EAs: The NSGA-II and Other Algorithms // IEEE Transactions on Evolutionary Computation. 2003. Vol. 7, no. 5. P. 503-515.

15Fortin, F-A., Grenier, S., Parizeau, M. Generalizing the Improved Run-time Complexity Algorithm for Non-dominated Sorting // Proceedings of Genetic and Evolutionary Computation Conference. 2013. P. 615-622; Buzdalov, M., Shalyto, A. A Provably Asymptotically Fast Version of the Generalized Jensen Algorithm for Non-Dominated Sorting//Parallel Problem Solving from Nature-PPSN XIII. 2014. P. 528537. (Lecture Notes in Computer Science ; 8672); Markina, M., Buzdalov, M. Towards Large-Scale Multiobjective Optimisation with a Hybrid Algorithm for Non-dominated Sorting // Parallel Problem Solving from Nature - PPSN XV, Vol. 1. 2018. P. 347-358. (Lecture Notes in Computer Science ; 11101).

боте Небро и Дурильо16, где было показано, что инкрементальная версия алгоритма NSGA-II имеет преимущества в качестве получаемых решений по сравнению с классической (при том же ограничении на число вычислений приспособленности особей), но существенно проигрывает по времени. Коллективом ученых под руководством Деба был предложен первый специализированный алгоритм недоминирующей сортировки, предназначенный для эффективного пересчета рангов при вставке новых точек, известный под названием ENLU (англ. Efficient Non-dominated Level Update)17, и хотя он все еще требует 0(п2к) времени в худшем случае на одну операцию вставки, на практике вставка чаще всего выполняется быстрее. Однако, несмотря на то, что данный результат был улучшен в последующей работе18, эффективность процедуры по добавлению точек и обновлению значений рангов в недоминирующей сортировке — инкрементальной недоминирующей сортировки — остается низкой. Кроме того, в указанных работах не рассматривался вопрос работы в контексте асинхронного вычисления приспособленности, когда множество потоков осуществляет вставку точек параллельно и независимо друг от друга.

Целью данной работы является повышение производительности инкрементальных многокритериальных эволюционных алгоритмов, в том числе при условии асинхронного вычисления приспособленности особей.

Для достижения указанной определены следующие задачи исследования:

а) уточнение задачи инкрементальной недоминирующей сортировки с целью повышения эффективности ее решения при условии сохранения совместимости с инкрементальными многокритериальными эволюционными алгоритмами;

б) разработка и экспериментальное исследование эффективных однопо-точных алгоритмов и структур данных для решения задач инкрементальной недоминирующей сортировки;

в) разработка и экспериментальное исследование многопоточных версий указанных алгоритмов, выработка рекомендаций по их использованию.

Объект исследования — процедура недоминирующей сортировки в инкрементальных многокритериальных эволюционных алгоритмах.

Предмет исследования — эффективные алгоритмы для инкрементальной недоминирующей сортировки, в том числе в контексте асинхронного вычисления функции приспособленности.

Соответствие паспорту специальности. Работа соответствует пункту 8 специальности «Модели и методы создания программ и программных систем

16Nebro, A. J., Durillo, J. J. On the Effect of Applying a Steady-State Selection Scheme in the Multi-Objective Genetic Algorithm NSGA-II // Nature-Inspired Algorithms for Optimisation. 2009. P. 435-456. (Studies in Computational Intelligence ; 193).

11 Li, K. [et al.]. Efficient Nondomination Level Update Method for Steady-State Evolutionary Multiob-jective Optimization // IEEE Transactions on Cybernetics. 2011. Vol. 41, no. 9. P. 2838-2849.

18Mishra, S., Mondal, S., Saha, S. Improved solution to the non-domination level update problem // Applied Soft Computing. 2011. Vol. 60. P. 336-362.

для параллельной и распределенной обработки данных, языки и инструментальные средства параллельного программирования» в части разработки эффективных алгоритмов асинхронной недоминирующей сортировки, предназначенных для создания распределенных многокритериальных эволюционных алгоритмов.

Основные положения, выносимые на защиту:

а) разработан алгоритм инкрементальной недоминирующей сортировки для случая к = 2, осуществляющий вставку одной точки за время O(n) в худшем случае и за время O(log n) в случае константного числа различных рангов точек, а также сопутствующая ему структура данных;

б) для случая к > 2 разработан алгоритм инкрементальной недоминирующей сортировки, основанный на использовании существующего алгоритма обычной недоминирующей сортировки, осуществляющий вставку одной точки за время O(n(logn)k-2);

в) разработаны два алгоритма для асинхронной недоминирующей сортировки, основанных на алгоритме инкрементальной недоминирующей сортировки — с использованием техники compare-and-set, а также с использованием блокировок, привязанных к рангам точек.

Научная новизна заключается в следующем:

а) в алгоритме инкрементальной недоминирующей сортировке для случая к = 2 используется новое, ранее не выявленное свойство, состоящее в том, что при вставке точки ранги остальных точек меняются следующим образом: среди всех точек, имеющих заданный ранг, ранг увеличивается у тех и только у тех точек, абсциссы которых попадают в некоторый отрезок;

б) в алгоритме инкрементальной недоминирующей сортировки для случая к > 2 обновление рангов точек производится итеративно от меньших рангов к большим, при этом на каждой итерации обновляются точки одного ранга, что позволило использовать свойство алгоритма обычной недоминирующей сортировки, согласно которому при константном числе рангов его время работы составляет O(n(logn)k-2);

в) алгоритмы для асинхронной недоминирующей сортировки ранее не предлагались, кроме того, в предложенных алгоритмах используются свойства инкрементальной недоминирующей сортировки для повышения их эффективности.

Методология и методы исследования: используются методы дискретной математики и вычислительной геометрии, объектно-ориентированное программирование, методы проведения и анализа вычислительных экспериментов.

Достоверность научных положений, выводов и доказательств, полученных в диссертации, подтверждается корректностью постановок задач, формальным доказательством утверждений относительно корректности и эффективности разработанных алгоритмов, результатами экспериментальных исследований.

Теоретическая значимость работы заключается в том, что в ней предложены алгоритмы инкрементальной недоминирующей сортировки, имеющие

на настоящий момент наилучшую асимптотическую сложность по сравнению с существующими алгоритмами, как потенциально применимыми, так и разработанными специально для данной задачи.

Практическая значимость работы заключается в существенном повышении эффективности инкрементальных многокритериальных эволюционных алгоритмов, использующих недоминирующую сортировку, а также их асинхронных вариантов.

Внедрение результатов работы. Результаты работы использованы при выполнении гранта Российского научного фонда (проект 17-71-20178 «Методы построения эффективных эволюционных алгоритмов», 2017-2020 гг.), а также в рамках государственной финансовой поддержки ведущих университетов Российской Федерации, субсидия 08-08 (НИР «Методы, модели и технологии искусственного интеллекта в биоинформатике, социальных медиа, киберфизических, биометрических и речевых системах», 2017-2019 гг.).

Результаты работы внедрены в программный продукт jMetal19, разработанный и поддерживаемый коллективом университета Малаги, г. Малага, Испания, что подтверждается письмом главного разработчика данного программного продукта А. Небро. Результаты работы также внедрены в учебный процесс факультета информационных технологий и программирования Университета ИТМО в курсе «Генетические и эволюционные вычисления» магистерской программы «Технологии разработки программного обеспечения», что подтверждается актом об использовании.

Апробация результатов работы. Основные результаты работы докладывались на следующих конференциях и семинарах:

а) Congress on Evolutionary Computation (2015, Сендай, Япония);

б) Genetic and Evolutionary Computation Conference (2015, Мадрид, Испания; 2017, Берлин, Германия; 2018, Киото, Япония);

в) Всероссийский конгресс молодых ученых (2015, 2018, Университет ИТМО).

Личный вклад автора. Автором выполнена разработка и реализация всех инкрементальных и асинхронных алгоритмов недоминирующей сортировки, выполнены все экспериментальные исследования. В работах, выполненных в соавторстве, Буздаловым М.В. осуществлялось общее руководство исследованием, кроме того, ему принадлежит реализация эффективного алгоритма обычной (не инкрементальной) недоминирующей сортировки, используемого в качестве составной части алгоритма инкрементальной недоминирующей сортировки для к > 2. Станкевичу А.С. в совместной с ним работе принадлежит первоначальная реализация инкрементального алгоритма NSGA-II.

Публикации. Основные результаты по теме диссертации изложены в четырех работах, опубликованных в международных изданиях, индексируемых в базе Scopus.

19https://github.com/jMetal/jMetal

Содержание работы

Во введении обосновывается актуальность исследований, проводимых в рамках данной диссертационной работы, приводится краткий обзор научной литературы по изучаемой проблеме, формулируется цель, ставятся задачи работы, излагается её научная новизна, теоретическая и практическая значимость.

Первая глава представляет собой аналитический обзор предметной области. Описываются различные постановки задач многокритериальной оптимизации и связанные с ними определения, такие как оптимальность по Парето, доминирование по Парето, множество Парето и фронт Парето. Мотивируется использование в диссертации постановки задачи многокритериальной оптимизации, требующей нахождения множества оптимальных по Парето решений, векторы значений критериев которых распределены равномерно по фронту Парето.

Для решения задач многокритериальной оптимизации в данной постановке часто применяются эволюционные алгоритмы. Причиной этому является то, что эволюционные алгоритмы в большинстве случаев работают с множествами решений-кандидатов (или, в терминах эволюционных алгоритмов, с популяциями особей), а следовательно, нахождение множества решений задачи многокритериальной оптимизации выполняется натуральным образом. В диссертации приведен краткий экскурс в общую историю эволюционных алгоритмов, а также более подробно приведен обзор наиболее важных классов многокритериальных эволюционных алгоритмов. Далее в содержании работы будет в основном использоваться терминология эволюционных алгоритмов: решение-кандидат задачи оптимизации называется особью, вектор значений критериев, соответствующих этому решению — приспособленностью данной особи, каждый критерий называется функцией приспособленности, рассматриваемое в настоящий момент множество решений называется популяцией.

Большинство современных многокритериальных эволюционных алгоритмов относится к одному из трех классов. Первый класс включает в себя алгоритмы, напрямую использующие отношение доминирования по Парето как основной инструмент ранжирования особей. Второй класс включает в себя алгоритмы, оптимизирующие некоторую функцию, принимающую на вход популяцию (или несколько популяций) и возвращающую число, характеризующее их качество; такие функции называются индикаторами. Третий класс алгоритмов сводит многокритериальную задачу оптимизации к серии однокритериальных подзадач оптимизации, которые решаются данными алгоритмами одновременно с возможностью миграции особей между решаемыми подзадачами и взаимного влияния соответствующих подзадачам популяций. Существуют также алгоритмы, сочетающие в себе черты нескольких классов. Отметим также, что ряд менее популярных, но тем не менее важных разновидностей эволюционных алгоритмов многокритериальной оптимизации — таких как, например, алгоритмы, интегрирующие предпочтения пользователей, а также алгоритмы интерактивной

оптимизации — нередко не могут быть отнесены ни к одному из этих классов, но в данной диссертации такие алгоритмы не рассматриваются.

Отношение доминирования по Парето используется алгоритмами из первого класса различным образом. Наиболее тривиальным из этих способов является поддержка архива особей, не доминируемых по Парето (далее просто «недоминируемых особей»). Также используется подсчет числа особей, доминирующих данную или доминируемых данной, и другие более сложные концепции.

Так, в алгоритмах SPEA20 и SPEA221 у каждой недоминируемой особи вычисляется сила Парето, зависящая от числа особей, которые доминируются данной особью, а для каждой из остальных особей вычисляется сумма значений силы Парето недоминируемых особей, доминирующих данную. Наконец, достаточно плодотворной идеей оказалось разделение популяции на «слои доминирования», известное как недоминирующая сортировка. Недоминирующая сортировка назначает каждой особи ранг по следующему правилу: всем недоминируемым особям в популяции назначается ранг 0, всем недоминируемым особям при условии удаления особей ранга 0 назначается ранг 1, всем недоминируемым особям при условии удаления особей ранга, меньшего г, назначается ранг г. Недоминирующая сортировка предложена в составе алгоритма NSGA22 и используется в ряде крайне популярных алгоритмов данного семейства (NSGA-II23, NSGA-Ш24), а также в значительном числе других алгоритмов.

Вычислительная сложность недоминирующей сортировки достаточно высока. В составе алгоритма NSGA для недоминирующей сортировки был предложен алгоритм со временной сложностью ©(п3к), где п — число особей, а к — число критериев оптимизации. В рамках алгоритма NSGA-П был предложен улучшенный алгоритм со сложностью <Э(п2к) по времени и 0(п2) по памяти. Так как квадратичная сложность все еще достаточно высока — в частности, это ограничивает размеры популяции, для которых недоминирующая сортировка может выполняться эффективно — в ряде последующих работ было разработано множество более эффективных алгоритмов. Большинство из них все еще имеет сложность <Э(п2к) по времени в худшем случае, однако на практике они работают гораздо быстрее. По состоянию на 2019 год, одними из

20Zitzler, E., Thiele, L. Multiobjective Evolutionary Algorithms: A Comparative Case Study and the Strength Pareto Approach // IEEE Transactions on Evolutionary Computation. 1999. Vol. 3, no. 4. P. 257271.

21Zitzler, E., Laumanns, M., Thiele, L. SPEA2: Improving the Strength Pareto Evolutionary Algorithm for Multiobjective Optimization // Proceedings of the EUROGEN'2001 Conference. 2001. P. 95-100.

22Srinivas, N., Deb, K. Multiobjective Optimization Using Nondominated Sorting in Genetic Algorithms // Evolutionary Computation. 1994. Vol. 2, no. 3. P. 221-248.

23Deb, K. [etal.]. A Fast and Elitist Multi-Objective Genetic Algorithm: NSGA-II//IEEE Transactions on Evolutionary Computation. 2002. Vol. 6, no. 2. P. 182-197.

24Deb, K., Jain, H. An Evolutionary Many-Objective Optimization Algorithm Using Reference-Point-Based Nondominated Sorting Approach, Part I: Solving Problems With Box Constraints // IEEE Transactions on Evolutionary Computation. 2013. Vol. 18, no. 4. P. 577-601.

лучших алгоритмов данного класса являются алгоритмы Best Order Sort25 и ENS-NDT26. Отдельно стоит также отметить такой алгоритм, как Merge Non-Dominated Sorting27, который демонстрирует эффективность лишь для больших значений числа критериев к. Существует также семейство алгоритмов, основанных на технике «разделяй-и-властвуй» со временем работы в худшем случае, составляющим O(n ■ (logn)k-1). Данное направление инициировано в работе Йен-сена28, затем последовал ряд теоретических и практических улучшений в работах Фортена с соавторами29 и Буздалова с соавторами30.

Для полноты картины также отметим, что недоминирующая сортировка не является самой вычислительно сложной процедурой, используемой в многокритериальных эволюционных алгоритмах — например, вычисление гиперобъема, популярного индикатора для оценки качества популяций, является #Р-трудной

Список литературы диссертационного исследования кандидат наук Якупов Илья Юрьевич, 2019 год

r □

8 □

A A

A A

Q

— y

ii i II 1 I o Q .......

1 10 100 Number of threads

Figure 33 - Five-dimensional DTLZ1 with artificial generation and sorting of 10000 floating-point values. Fine grained locking algorithm is marked with blue squares, CAS-based algorithm - with red triangles, ENLU - with yellow circles, coarse-grained - with black circles. Generational NSGA-II (green circles) and single-threaded incremental algorithm (magenta triangles) are also included in the

comparison.

able problem sizes. It can possibly be done by checking in O(nk) whether the regions dominated by two different sets of moving points, that are manipulated by different threads, intersect in the current level: if they do not, then this level, and all subsequent levels, can be processed by these two threads only with minor resource sharing, since the most expensive parts will operate with non-intersecting sets of points.

Yet another possibility, which may find its use in heterogeneous computing systems (where the internals of an evolutionary algorithm is run with different computation resources than fitness evaluation) is a special flavor of an asynchronous algorithm which, on the arrival of a fitness thread, hands the next task immediately, and only then inserts the evaluated point into the data structure. This should reduce the idle rate of fitness-related computation resources, which are typically more expensive. However, an impact of this design on the convergence of an algorithm, as compared to the one implemented in this paper, may be non-trivial and should be carefully investigated.

The presented results have been published in [47].

Conclusion

The main results of the work are as follows.

- An algorithm for incremental non-dominated sorting is proposed for the k = 2 case, which performs insertion of a single point in the O(n) worst-case time and in O (log n) time when the number of different ranks of points is constant. It uses a novel property, which constrains the possible rank changes during insertion of a single point in the following way: among all points that have a given rank, those and only those points have their rank increased whose abscissa values fall into a certain segment. Experimental evaluation shows that this algorithm performs noticeably better than the existing solutions to this problem.

- An algorithm for incremental non-dominated sorting is proposed for the k > 2 case, which performs insertion of a single point in the O(n(log n)k-2) time, which internally uses an existing algorithm for the conventional non-dominated sorting. It updates point ranks iteratively from smaller to greater ranks, while on each iteration only points with the same rank are updated, which enables using the following property of the existing algorithm for the conventional non-dominated sorting: when the number of ranks is constant, its running time is O(n(log n)k-2). Experimental evaluation suggests that this algorithm typically performs better than the existing solutions, which was supported by statistical testing.

- Two algorithms for asynchronous incremental non-dominated sorting are proposed, which are based on the single-threaded algorithm, where the first one employs the compare-and-set technique, and the second one uses fine-grained locks attached to the point sets having the same rank. Such algorithms have never been proposed before, and certain properties of the solved problem are extensively used in the proposed algorithms to enhance their efficiency. Experimental evaluation suggests that the algorithm using fine-grained locking is the most efficient for practically reasonable numbers of threads.

References

1. Luke, S. Essentials of Metaheuristics. — Lulu, 2009. — 253 p.

2. Deb, K. [et al.]. A Fast and Elitist Multi-Objective Genetic Algorithm: NSGA-II // IEEE Transactions on Evolutionary Computation. — 2002. — Vol. 6, no. 2.—P. 182-197.

3. Zitzler, E., Laumanns, M., Thiele, L. SPEA2: Improving the Strength Pareto Evolutionary Algorithm for Multiobjective Optimization // Proceedings of the EU-R0GEN'2001 Conference. —2001. — P. 95-100.

4. Corne, D. W. [et al.]. PESA-II: Region-based Selection in Evolutionary Multi-objective Optimization//Proceedings of Genetic and Evolutionary Computation Conference. — 2001. — P. 283-290.

5. Zitzler, E., Kunzli, S. Indicator-Based Selection in Multiobjective Search // Parallel Problem Solving from Nature - PPSN VIII. — 2004. — P. 832-842. — (Lecture Notes in Computer Science ; 3242).

6. Bader, J., Zitzler, E. HypE: An Algorithm for Fast Hypervolume-Based Many-Objective Optimization // Evolutionary Computation. — 2011. — Vol. 19, no. 1.—P. 45-76.

7. Zhang, Q., Li, H. MOEA/D: A Multiobjective Evolutionary Algorithm Based on Decomposition // IEEE Transactions on Evolutionary Computation. — 2007. — Vol. 11, no. 6.—P. 712-731.

8. Deb, K., Jain, H. An Evolutionary Many-Objective Optimization Algorithm Using Reference-Point-Based Nondominated Sorting Approach, Part I: Solving Problems With Box Constraints // IEEE Transactions on Evolutionary Computation. — 2013. — Vol. 18, no. 4. — P. 577-601.

9. Nebro, A. J., Durillo, J. J. On the Effect of Applying a Steady-State Selection Scheme in the Multi-Objective Genetic Algorithm NSGA-II // Nature-Inspired Algorithms for Optimisation. — 2009. — P. 435-456. — (Studies in Computational Intelligence ; 193).

10. Srinivas, N., Deb, K. Multiobjective Optimization Using Nondominated Sorting in Genetic Algorithms // Evolutionary Computation. — 1994. — Vol. 2, no. 3. — P. 221-248.

11. Roy, P. C., Islam, M.M., Deb, K. Best Order Sort: A New Algorithm to Non-dominated Sorting for Evolutionary Multi-objective Optimization // Proceedings of Genetic and Evolutionary Computation Conference Companion. — 2016. — P. 1113-1120.

12. Roy, P. C., Deb, K., Islam, M.M. An Efficient Nondominated Sorting Algorithm for Large Number of Fronts // IEEE Transactions on Cybernetics. — 2019. — Vol. 49,no. 3.— P. 859-869.

13. Gustavsson, P., Syberfeldt, A. A New Algorithm Using the Non-dominated Tree to Improve Non-dominated Sorting // Evolutionary Computation. — 2018. — Vol. 26, no. 1.—P. 89-116.

14. Zhang, X.[et al.]. An Efficient Approach to Nondominated Sorting for Evolutionary Multiobjective Optimization//IEEE Transactions on Evolutionary Computation. — 2015. — Vol. 19, no. 2. — P. 201-213.

15. Jensen, M. T. Reducing the Run-time Complexity of Multiobjective EAs: The NSGA-II and Other Algorithms // IEEE Transactions on Evolutionary Computation. — 2003. — Vol. 7, no. 5. — P. 503-515.

16. Fortin, F.-A., Grenier, S., Parizeau, M.Generalizing the Improved Run-time Complexity Algorithm for Non-dominated Sorting // Proceedings of Genetic and Evolutionary Computation Conference. — 2013. — P. 615-622.

17. Buzdalov, M., Shalyto, A. A Provably Asymptotically Fast Version of the Generalized Jensen Algorithm for Non-Dominated Sorting // Parallel Problem Solving from Nature - PPSN XIII. — 2014. — P. 528-537. — (Lecture Notes in Computer Science ; 8672).

18. Markina, M., Buzdalov, M.Towards Large-Scale Multiobjective Optimisation with a Hybrid Algorithm for Non-dominated Sorting // Parallel Problem Solving from Nature - PPSN XV, Vol. 1. — 2018. — P. 347-358. — (Lecture Notes in Computer Science ; 11101).

19. Li, K. [et al.]. Efficient Nondomination Level Update Method for Steady-State Evolutionary Multiobjective Optimization // IEEE Transactions on Cybernetics. — 2017. — Vol. 47, no. 9. — P. 2838-2849.

20. Mishra, S., Mondal, S., Saha, S. Improved solution to the non-domination level update problem // Applied Soft Computing. — 2017. — Vol. 60. — P. 336-362.

21. Branke, J. [et al.]. Multiobjective Optimization - Interactive and Evolutionary Approaches. — Springer, 2008.

22. Wikipedia. Эффективность по Парето. — URL: https://ru.wikipedia.org/wiki/ %D0%AD%D1%84%D1%84%D0%B5%D0%BA%D1%82%D0%B8%D0% B2%D0%BD%D0%BE%D1 %81 %D1 %82%D1 %8C_%D0%BF%D0%BE_ %D0%9F%D0%B0%D1%80%D0%B5%D1%82%D0%BE.

23. Brockhoff, D., Wagner, T. GECCO 2016 Tutorial on Evolutionary Multiobjective Optimization // Proceedings of Genetic and Evolutionary Computation Conference Companion. — 2016. — P. 201-227.

24. Schaffer, D. Multiple Objective Optimization with Vector Evaluated Genetic Algorithms // Proceedings of the 1st International Conference on Genetic Algorithms. — 1985. — P. 93-100.

25. Fonseca, C. M., Fleming, P. J. Nonlinear System Identification with Multiobjective Genetic Algorithm // Proceedings of the World Congress of the International Federation of Automatic Control. — 1996. — P. 187-192.

26. Horn, J., Nafpliotis, N., D.E.Goldberg. A Niched Pareto Genetic Algorithm for Multiobjective Optimization//Proceedings of IEEE Conference on Evolutionary Computation. Vol. 1. — 1994. — P. 82-87.

27. Erickson, M., Mayer, A., Horn, J. The Niched Pareto Genetic Algorithm 2 Applied to the Design of Groundwater Remediation Systems // Proceedings of International Conference on Evolutionary Multi-Criterion Optimization. — 2001. — P. 681-695. — (Lecture Notes in Computer Science ; 1993).

28. Zitzler, E., Thiele, L. Multiobjective Evolutionary Algorithms: A Comparative Case Study and the Strength Pareto Approach // IEEE Transactions on Evolutionary Computation. — 1999. — Vol. 3, no. 4. — P. 257-271.

29. Li, K. [et al.]. Efficient Non-domination Level Update Approachfor Steady-State Evolutionary Multiobjective Optimization : tech. rep. / Michigan State University. — 2014. — COIN 2014014. — URL: www.egr.msu.edu/~kdeb/papers/ c2014014.pdf.

30. Yagoubi, M., Thobois, L., Schoenauer, M.Asynchronous evolutionary multi-objective algorithms with heterogeneous evaluation costs // Proceedings of IEEE Congress on Evolutionary Computation. —2011. —P. 21-28.

31. Wang, Q. [et al.]. Parameterization of NSGA-II for the Optimal Design of Water Distribution Systems // Water. — 2019. — Vol. 11(5):971.

32. Kung, H.-T., Luccio, F., Preparata, F. P. On Finding the Maxima of a Set of Vectors // Journal of ACM. — 1975. — Vol. 22, no. 4. — P. 469-476.

33. Churchill, A. W., Husbands, P, Philippides, A. Tool Sequence Optimization using Synchronous and Asynchronous Parallel Multi-Objective Evolutionary Algorithms with Heterogeneous Evaluations // Proceedings of IEEE Congress on Evolutionary Computation. — 2013. — P. 2924-2931.

34. Durillo, J. J.[et al.]. A study of master-slave approaches to parallelize NSGA-II // Proceedings of IEEE International Symposium on Parallel and Distributed Processing. — 2008. — P. 1-8.

35. Yagoubi, M., Schoenauer, M.Asynchronous master/slave MOEAs and heterogeneous evaluation costs // Proceedings of Genetic and Evolutionary Computation Conference. — 2012. — P. 1007-1014.

36. Yakupov, I., Buzdalov, M. Incremental Non-Dominated Sorting with O(N) Insertion for the Two-Dimensional Case //Proceedings of IEEE Congress on Evolutionary Computation. — 2015. — P. 1853-1860.

37. Vuillemin, J. A Unifying Look at Data Structures // Communications of ACM. — 1980. — Vol. 23, no. 4. — P. 229-239.

38. Sleator, D. D., Tarjan, R. E. Self-Adjusting Binary Search Trees // Journal of ACM. — 1985. — Vol. 32, no. 3. — P. 652-686.

39. Zitzler, E., Deb, K., Thiele, L. Comparison of Multiobjective Evolutionary Algorithms: Empirical Results // Evolutionary Computation. — 2000. — Vol. 8, no. 2.—P. 173-195.

40. Deb, K. [et al.]. Scalable Test Problems for Evolutionary Multiobjective Optimization //Evolutionary Multiobjective Optimization. Theoretical Advances and Applications. —2005. — P. 105-145.

41. Huband, S. [et al.]. A review of multiobjective test problems and a scalable test problem toolkit // IEEE Transactions on Evolutionary Computation. — 2006. — Vol. 10, no. 5.—P. 477-506.

42. Buzdalov, M., Yakupov, I., Stankevich, A. Fast Implementation of the Steady-State NSGA-II Algorithm for Two Dimensions Based on Incremental Non-Dominated Sorting // Proceedings of Genetic and Evolutionary Computation Conference. — 2015. — P. 647-654.

43. Yakupov, I., Buzdalov, M.Improved Incremental Non-dominated Sorting for Steady-State Evolutionary Multiobjective Optimization // Proceedings of Genetic and Evolutionary Computation Conference. — 2017. — P. 649-656.

44. Hoare, C. Monitors: an operating system structuring concept // Communications of ACM. — 1974. — Vol. 17, no. 10. — P. 549-557.

45. Kogan, A., Petrank, E. Wait-free queues with multiple enqueuers and dequeuers // Proceedings of the ACM SIGPLAN Symposium on Principles and Practice of Parallel Programming. — 2011. — P. 223-234.

46. Hadka, D. MOEA Framework - A Free and Open Source Java Framework for Multiobjective Optimization. Version 2.12. — 2015. — URL: http://www. moeaframework.org.

47. Yakupov, I., Buzdalov, M. On asynchronous non-dominated sorting for steady-state multiobjective evolutionary algorithms // Proceedings of Genetic and Evolutionary Computation Conference Companion. — 2018. — P. 205-206. — URL: https://arxiv.org/pdf/1804.05208.

List of Figures

Р.1 Вставка точки в множество двумерных отсортированных точек . 19

Р.2 Структура данных «дерево деревьев» для двумерной инкрементальной недоминирующей сортировки ............... 20

Р.3 Время работы вставки N точек, случайно сгенерированных в

квадрате ................................. 21

Р.4 Вставка одной точки с послойным обновлением рангов обычной

недоминирующей сортировкой .................... 24

Р.5 Пример результатов экспериментальных исследований: задача

DTLZ7, K = 3............................. 24

Р.6 Ускорение обработки вставки точек, достигаемое исследуемыми алгоритмами, в случае моделируемого времени вычисления приспособленности, равного 1 мс .................... 26

5.1 Insertion of a point in a set of 2D sorted points............. 39

5.2 The "tree-of-trees" data structure for two-dimensional incremental non-dominated sorting......................... 40

5.3 The wall-clock time for insertion of N points randomly sampled from

a square................................. 41

5.4 Insertion of one point with layer-wise rank updates by offline non-dominated sorting ............................ 43

5.5 An example of experimental results: the DTLZ7 problem, K = 3 . . 44

5.6 Speedup of point insertion, reachable by the proposed algorithms in

the case of modeled fitness evaluation time of 1 ms.......... 45

1 An example of a multiobjective problem with a Pareto front ..... 54

2 Non-dominated sorting: An example.................. 58

3 The original NSGA-II non-dominated sorting procedure. For each point the algorithm stores the set of dominated points, as well as the count of individuals that dominate over it. If this count reaches zero, the point should be assigned to the current non-domination level. The overall running time complexity is O(KN2).............. 61

4 Illustration of the SweepA procedure.................. 64

5 Example of an ENLU behavior. The darkest point is the newly inserted point. All other gray points are the points which change their rank. The gray lines indicate the comparisons which ENLU performs. The initial

levels are shown on the left, the resulting levels on the right...... 67

6 The main procedure of ENLU. At first, the rank of the added point is determined. After that, the point is being added to the corresponding non-domination level. Here, P is the sorted population, xc is the added individual, and l is the maximal known rank of any individual in the population. Fi is the i-th ranked non-domination level......... 68

7

8

9

10

11

12

13

14

15

16

17

18

19

20

21

22

23

24

The sub-procedure of ENLU that pushes the "evicted" points to the "furthef' levels. Here, S is the set of the "pushed" individuals, and i is the index of the level, to which the points of S should be added . .

A pseudocode for the data structure ..................

Data structure for the algorithm - the "tree of trees". Nodes of the "high-level" tree correspond to the levels. Each level is, in turn, represented by a "low-level" tree, where nodes are sorted by the first objective. In each node of the "high-level" tree, the subtree size is stored as well (the numbers in parentheses shown near to the level number). In each node of the "low-level" tree we store the associated random value, required by the Cartesian tree (the numbers in parentheses shown near to the level number). Note that level numbers are not stored in nodes

explicitly, they are just shown for convenience.............

A pseudocode for determining the smallest level which doesn't dominate the given solution .........................

A pseudocode for helper procedures in a high-level tree .......

A pseudocode for insertion of a solution into a high-level tree .... A pseudocode for querying the fc-th solution in lexicographical order . A pseudocode for querying the fc-th solution in lexicographical order,

continued ................................

The "square" test: Illustration and experimental results ........

The "parallel" test: Illustration and experimental results ........

The "diag1" test: Illustration and experimental results.........

The "diag2" test: Illustration and experimental results.........

The "parper" test: Illustration and experimental results ........

The working principles of incremental non-dominated sorting for K > 2 Running times of single point insertions, in microseconds, of the algorithms on two-dimensional ZDT1-ZDT4 and ZDT6. Abscissa marks are the number of generations simulated. The algorithms, left to right:

DC, DC1, DC2, Level, ENLU, 2D....................

Running times of single point insertions, in microseconds, of the algorithms on three-dimensional DTLZ1-DTLZ4 and DTLZ7. Abscissa marks are the number of generations simulated. The algorithms, left

to right: DC, DC1, DC2, Level, ENLU.................

Running times of single point insertions, in microseconds, of the algorithms on DTLZ1 with dimensions 3, 4, 6, 8, 10. Abscissa marks are the number of generations simulated. The algorithms, left to right:

DC, DC1, DC2, Level, ENLU......................

Running times of single point insertions, in microseconds, of the algorithms on DTLZ2 with dimensions 3, 4, 6, 8, 10. Abscissa marks are the number of generations simulated. The algorithms, left to right: DC, DC1, DC2, Level, ENLU......................

25 A lock-based implementation of the compare-and-set function to illustrate its semantics............................. 106

26 Wall-clock running times of population modification procedures of the steady-state NSGA-II algorithm. FNDS and INDS are single-threaded, all other algorithms use 3, 6 and 12 threads, depicted left

to right.................................. 113

27 Resulting hypervolumes of the steady-state NSGA-II algorithm. FNDS and INDS are single-threaded, all other algorithms use 3, 6 and

12 threads, depicted left to right..................... 115

28 The wall clock running time reduction (compared to a single-thread configuration), achieved by adding worker threads. The value of speedup is the average single-threaded running time, divided by the running time of the current multi-threaded execution. The evaluated function is the three-dimensional DTLZ1 with artificial delay of one millisecond. Red markers are the results, shown by the CAS-based algorithm, black markers correspond to the coarse-grained locking algorithm and the blue ones represent the fine-grained locking algorithm. 117

29 The wall clock running time reduction (compared to a single-thread configuration), achieved by adding worker threads. The value of speedup is the average single-threaded running time, divided by the running time of the current multi-threaded execution. The evaluated function is the three-dimensional DTLZ1 with artificial delay of five milliseconds. Red markers are the results, shown by the CAS-based algorithm, black markers correspond to the coarse-grained locking algorithm and the blue ones represent the fine-grained locking algorithm. 118

30 The wall clock running time reduction (compared to a single-thread configuration), achieved by adding worker threads. The value of speedup is the average single-threaded running time, divided by the running time of the current multi-threaded execution. The evaluated function is the four-dimensional DTLZ1 with artificial delay of one millisecond. Red markers are the results, shown by the CAS-based algorithm, black markers correspond to the coarse-grained locking algorithm and the blue ones represent the fine-grained locking algorithm. 119

31 Three-dimensional DTLZ1 with artificial generation and sorting of 10000 floating-point values. Fine grained locking algorithm is marked with blue squares, CAS-based algorithm - with red triangles, ENLU -with yellow circles, coarse-grained - with black circles. Generational NSGA-II (green circles) and single-threaded incremental algorithm (magenta triangles) are also included in the comparison......... 120

32 Four-dimensional DTLZ1 with artificial generation and sorting of 10000 floating-point values. Fine grained locking algorithm is marked with blue squares, CAS-based algorithm - with red triangles, ENLU -with yellow circles, coarse-grained - with black circles. Generational NSGA-II (green circles) and single-threaded incremental algorithm (magenta triangles) are also included in the comparison......... 121

33 Five-dimensional DTLZ1 with artificial generation and sorting of 10000 floating-point values. Fine grained locking algorithm is marked with blue squares, CAS-based algorithm - with red triangles, ENLU -with yellow circles, coarse-grained - with black circles. Generational NSGA-II (green circles) and single-threaded incremental algorithm (magenta triangles) are also included in the comparison......... 122

List of Tables

Р.1 Результаты экспериментов в двумерном случае, задача ZDT1 ... 22

S.1 Experiment results for two dimensions, problem ZDT1........ 42

1 Results of experiments. Notation: "time" stands for the running time (in seconds), "cmp" for the number of comparisons, "med" for the median, "IQR" for the interquartile range. "INDS" is the proposed implementation, "ENLU" is the ENLU-based implementation, "debNDS" is the standard implementation which uses the Deb's fast non-dominated sorting, "gen" corresponds to the generational variant of the NSGA-II algorithm and "ss" to the steady-state variant. Gray cells denote minima of the running time and of the number of comparisons (for generational and steady-state variants separately)........... 90

2 Settings for population sizes and computational budgets........ 112

Appendix A Copies of Author's Publications

1. Yakupov, I., Buzdalov, M. On asynchronous non-dominated sorting for steady-state multiobjective evolutionary algorithms // Proceedings of Genetic and Evolutionary Computation Conference Companion. — 2018. — P. 205-206. — URL: https://arxiv.org/pdf/1804.05208.

2. Yakupov, I., Buzdalov, M.Improved Incremental Non-dominated Sorting for Steady-State Evolutionary Multiobjective Optimization // Proceedings of Genetic and Evolutionary Computation Conference. — 2017. — P. 649-656.

3. Buzdalov, M., Yakupov, I., Stankevich, A. Fast Implementation of the Steady-State NSGA-II Algorithm for Two Dimensions Based on Incremental Non-Dominated Sorting // Proceedings of Genetic and Evolutionary Computation Conference. — 2015. — P. 647-654.

4. Yakupov, I., Buzdalov, M.Incremental Non-Dominated Sorting with O(N) Insertion for the Two-Dimensional Case // Proceedings of IEEE Congress on Evolutionary Computation. — 2015. — P. 1853-1860.

Fast Implementation of the Steady-State NSGA-II Algorithm for Two Dimensions Based on Incremental Non-Dominated Sorting

Maxim Buzdalov

ITMO University 49 Kronverkskiy ave. Saint-Petersburg, Russia

[e-mail]

Ilya Yakupov

ITMO University 49 Kronverkskiy ave. Saint-Petersburg, Russia

[e-mail]

Andrey Stankevich

ITMO University 49 Kronverkskiy ave. Saint-Petersburg, Russia

stankev@g mail.com

ABSTRACT

Genetic algorithms (GAs) are widely used in multi-objective optimization for solving complex problems. There are two distinct approaches for GA design: generational and steady-state algorithms. Most of the current state-of-the-art GAs are generational, although there is an increasing interest to steady-state algorithms as well. However, for algorithms based on non-dominated sorting, most of steady-state implementations have higher computation complexity than their generational counterparts, which limits their applicability.

We present a fast implementation of a steady-state version of the NSGA-II algorithm for two dimensions. This implementation is based on a data structure which has O(N) complexity for single solution insertion and deletion in the worst case. The experimental results show that our implementation works noticeably faster than steady-state NSGA-II implementations which use fast non-dominated sorting.

CCS Concepts

•Theory of computation ^ Mathematical optimization; Data structures design and analysis;

Keywords

steady-state, multi-objective, NSGA-II, non-dominated sorting, incremental updates

1. INTRODUCTION

In the K-dimensional space, a point A = (ai,...,ax) is said to dominate a point B = (bi,..., bx) when for all 1 < i < K it holds that < bi and there exists j such that aj < bj. Non-dominated sorting of points in the K-dimensional space is a procedure of marking all points which are not dominated by any other point with the rank of 0, all points which are dominated by at least one point of the rank of 0 are marked with the rank of 1, all points which are

Permission to make digital or hard copies of part or all of this work for personal or classroom use is granted without fee provided that copies are not made or distributed for profit or commercial advantage and that copies bear this notice and the full citation on the first page. Copyrights for third-party components of this work must be honored. For all other uses, contact the owner/author(s). GECCO '15 July 11-15, 2015, Madrid, Spain © 2015 Copyright held by the owner/author(s). ACM ISBN 978-1-4503-3472-3/15/07. DOI: http://dx.doi.org/10.1145/2739480.2754728

dominated by at least one point of the rank i — 1 are marked with the rank of i.

Many well-known and widely used multi-objective evolutionary algorithms use the procedure of non-dominated sorting, or the procedure of determining the non-dominated solutions, which can be reduced to non-dominated sorting. These algorithms include NSGA-II [5], PESA [4], PESA-II [3], SPEA2 [18], PAES [10], PDE [1], and many more. The time complexity of a single iteration of these algorithms is often dominated by the complexity of a non-dominated sorting algorithm, so optimization of the latter makes such multi-objective evolutionary algorithms faster.

In Kung et al [11], the algorithm for determining the non-dominated solutions is proposed with the running time complexity of O(N logx-i N), where N is the number of points and K is the dimension of the space. It is possible to use this algorithm to perform non-dominated sorting: first, the non-dominated solutions are found and assigned the rank of 0. Then, these solutions are removed, the non-dominated solutions from the remaining ones are found and assigned the rank of 1. The process repeats until all the solutions are removed. This leads to the complexity of O(N2 logx-i N) in the worst case, if the maximum rank of a point in the result is O(N).

Jensen [9] was the first to propose an algorithm for non-dominated sorting with the complexity of O(Nlogx-i N). However, his algorithm was developed for the assumption that no two points share a common value for any objective, and the complexity was proven for the same assumption. The first attempt to fix this issue belongs, to the best of the authors' knowledge, to Fortin et al [7]. The corrected (or, as in [7], "generalized") algorithm works in all cases, and for the general case the performance is still O(Nlogx-i N), but the only upper bound that was proven for the worst case is O(N2K). Finally, Buzdalov et al in [2] proposed several modifications to the algorithm of Fortin et al to make the O(Nlogx-i N) bound provable as well.

Evolutionary algorithms have a big advantage due to their great degree of parallelism, however, synchronous variants (which wait for evaluation of all individuals, then recompute their internal state) have only a limited applicability for distributed systems. Even on multicore computers an algorithm may have a poor performance if it spends big periods of time between fitness evaluations without using most of computer resources. To overcome these limitations, steady-state algorithms are developed, often with an intention to

become asynchronous. Particularly, a steady-state version of the NSGA-II algorithm was developed [13] which showed good convergence rate and high quality of Pareto front approximation. However, the running time of it is poor.

It is possible to perform incremental non-dominated sorting by doing a complete non-dominated sorting from scratch every time an element is added. However, running times become very high: O(KN2) for a single insertion when the fast non-dominated sorting [5] is used, or O(N logK-1 N) when the sorting from [2] is used. Thus, it is needed to develop new algorithms and data structures to handle incremental non-dominated sorting efficiently.

However, almost no such algorithms and data structures have been developed so far. To our best knowledge, the only paper which addresses this issue is a technical report by Li et al. [12]. In that report, a procedure called "Efficient Non-domination Level Update" is introduced, which has the complexity of O(NKVN) for a single insertion when solutions are spread evenly over layers. This procedure was shown experimentally to be quite efficient, however, the worst-case complexity for a single insertion is still O(N2K).

This paper briefly describes the incremental non-dominated sorting algorithm for the two-dimensional case (K = 2) which has worst-case insertion and deletion complexity of O(N) and worst-case query complexity of O(log N). It is described in our previous paper [16] in more detail. We additionally modified it to support evaluation of the crowding distance, a density measure employed by the NSGA-II algorithm. After that, the steady-state version of the NSGA-II algorithm is presented. Various implementation details are discussed, and results of experiments are presented.

2. INCREMENTAL NON-DOMINATED SORTING ALGORITHM

In this section, we briefly describe the proposed incremental non-dominated sorting algorithm and the data structure which supports it. They are described in detail in [16] along with necessary proofs.

The choice of the Cartesian tree as an underlying balanced search tree is discussed in Section 2.1. The data structure design is described in Section 2.2. The procedure of solution lookup (finding which layer a solution belongs to) is described in Section 2.3. The procedure of solution insertion is described in Section 2.4. The procedure of querying the k-th solution (in some predefined order) together with its layer number and crowding distance is described in Section 2.5. The worst solution deletion is described in Section 2.6.

When discussing the runtime analysis, we denote by N the total number of solutions stored in the data structure and by M the current number of non-domination layers. For the sake of brevity, non-domination layers are called just "layers" in the rest of the paper.

2.1 Cartesian Trees

To implement the algorithm, we need to have a data structure for container of elements which performs the following operations in O(log N):

• search of an element in the container;

• split of the container by key into two parts (the elements less than the key and the elements not less than the key);

1: structure Solution

2: - a solution to the optimization problem 3: X : Real - the first objective 4: Y : Real - the second objective 5: end structure 6: structure LLNode 7: - a node of a low-level tree 8: L : LLNode - the left child 9: R : LLNode - the right child 10: V : Solution - the node key 11: P : LLNode - the previous-in-order node 12: N : LLNode - the next-in-order node 13: S : Integer - the subtree size 14: end structure 15: structure HLNode 16: - a node of a high-level tree 17: L : HLNode - the left child 18: R : HLNode - the right child 19: P : HLNode - the previous-in-order node 20: N : HLNode - the next-in-order node 21: V : LLNode - the node key 22: S : Integer - the subtree size 23: W : Integer - the sum of key sizes in the subtree 24: end structure

Figure 1: A pseudocode for the data structure

• merge of two containers C1 and C2 (every element from C1 is not greater than every element from C2).

We will call data structures which fulfil these requirements "split-merge balanced search trees". There are several such data structures, including Cartesian Tree [15] and Splay Tree [14]. In the case of Cartesian Tree, the O(log N) bound holds with high probability, while Splay Tree has amortized O(log N) bounds. From the mentioned data structures, Cartesian Tree generally performs slightly better in practice, so we use it in an implementation of our algorithm.

2.2 Data Structure

The idea of the data structure is to arrange layers in a binary search tree (each tree node corresponds to a layer) in the increasing order of their numbers. Each layer, in turn, is represented by a binary search tree itself, where solutions are sorted in the increasing order of their first objective. Since for two different solutions a and b from the same layer it holds that either ax > bx and ay < by or ax < bx and ay > by, solutions in each layer are effectively sorted in the decreasing order of their second objective as well.

For the sake of brevity, we denote the tree of layers as the "high-level" tree and every tree containing layer elements as a "low-level" tree. The pseudocode for the data structure is given in Fig. 1. The resulting data structure is presented in Fig. 2.

Note that the high-level tree can be an ordinary balanced tree, while every low-level tree should be a split-merge tree. However, to evaluate the number of a certain layer in O(log M), one needs to store the number of tree elements in a subtree in each node of the high-level tree. Additionally, to move between adjacent layers in O(1), nodes should be augmented with pointers to the previous-in-order and the next-in-order nodes (which can be done without affecting O(log N) performance of basic operations).

Figure 2: The data structure for the algorithm — the "tree of trees". Nodes of the "high-level" tree correspond to the layers. Each layer is, in turn, represented by a "low-level" tree, where nodes are sorted by the first objective. Note that layer numbers are not stored in nodes explicitly, they are just shown for convenience.

2.3 Lookup

Given a low-level tree T and a solution s, it is possible to find if s is dominated by at least one solution from T in O(log |T|). To do it, one needs to find a solution u from T such that ux < sx and ux is maximum possible, which can be done by traversing the tree T from its root. If u is found and dominates s, then a dominating solution from T is found, otherwise, no solution from T dominates s.

Using this algorithm, one can traverse the high-level tree and find a layer with the minimum number which does not dominate a certain solution s. The algorithm is presented in Fig. 3.

A rough estimation of the running time is O(log M log N), where O(log M) is an estimation of the height of the highlevel tree, and O(log N) is an estimation of heights of all low-level trees.

However, one can perform a better estimation using the following idea. There are k = O(log M) layers which were tested for domination. Let their sizes be Li... Lk, and Li + ... + Lk < N. The running time for a layer of size Li can be expressed as O(1 + log Li) (we add extra 1 to handle a condition of log Li = o(1)). The total running time of a single lookup is:

O k + J]log Li

Due to Cauchy's inequality, ^i=i log Li < k log(N/k), which finally gives the following complexity of a lookup operation:

O

g M 1+1

N

which, due to the fact that M < N and log(N/ log M) is w(1), can be simplified to:

O ( log M log -

gM

When N is fixed and M varies, this expression reaches its maximum at M = ©(N), yielding O((logN)2) worst-case running time.

8 B — T

9 T i— T.R

10 else

11 T — T.L

12 : end if

function LowLevelDominates(T, s)

- returns whether any solution from T dominates s T : LLNode - the root node of the low-level tree s : Solution - the solution to test for domination B — null - the best node so far while T = null do if T.V.X < s.X then

end while if B = null then return false end if

return B.Y < s.Y or B.Y = s.Y and B.X < s.X end function

function SmallestNonDominatingLayer(H, s)

- returns the layer with the smallest index from H

- which does not dominate s

H : HLNode - the root node of the high-level tree s : Solution - the solution to find a layer for I — 0 - the number of dominating layers so far B — null - the best node so far while H = null do

if LowLevelDominates(H.V, s) then I — I + H.S H — H.R if H = null then

I — I - H.S end if else

B—H H — H.L end if end while return (B, I) end function

Figure 3: A pseudocode for determining the smallest layer which doesn't dominate the given solution

M

N

Figure 4: An example of the insertion process. Solutions which don't change their layer nodes (not numbers!) during the insertion process are white. A solution which is being inserted is black. Two clusters of solutions which together change their layer node are dark-gray and light-gray, correspondingly.

2.4 Insertion

Given a high-level tree H and a solution s, the insertion procedure updates H so that s is included in one of its low-level trees.

A key idea of fast implementation of insertion procedure is the fact that solutions who change their layers form contiguous pieces in their original layers and remain contiguous in their new layers as well. Fig. 4 illustrates an example of the insertion process.

The algorithm for a solution insertion is given in Fig. 5. It maintains a low-level tree which, at each stage, contains the solutions which needs to be inserted to the next layer. Initially it consists of the single solution which needs to be inserted. The layer to insert is initially found using the performing the SmallestNonDominatingLayer operation from Fig. 3.

The insertion algorithm works in iterations, each iteration pushes solutions to the layer that is immediately dominated by the layer of the previous iteration. On each iteration, the following operations are performed:

• The low-level tree of the current layer is split in three parts using the current pushed set of solutions C in the following way:

— the "left part" TL consists of all solutions from the current layer whose X coordinates are less than the smallest X coordinate of a solution from C;

— the "middle part" Tm consists of all solutions from the current layer which are dominated by at least one solution from C;

— the "right part" Tr consists of all solutions from the current layer whose Y coordinates are less than the smallest Y coordinate of a solution from C.

• The current layer is built by merging the trees TL, C and Tr.

• If both Tl and Tr are empty, this means that the entire level was dominated by solutions from C. In turn, this means that a new layer consisting entirely of Tm should be inserted just after the current level. All remaining layers will effectively have their index increased by one. The insertion procedure stops here.

1: function SplitX(T, s) 2: - splits a tree T into two trees L, R 3: - such that for all l e L holds l.X < s.X 4: - and for all r e R holds r.X > s.X 5: T : LLNODE 6: s : Solution 7: end function 8: function SplitY(T, s) 9: - splits a tree T into two trees L, R 10: - such that for all l e L holds l.Y > s.Y 11: - and for all r e R holds r.Y < s.Y 12: T : LLNODE 13: s : Solution 14: end function 15: function Merge(L, R)

16: - merges two trees L and R into a single one

17: - given for any l e L and r e R holds l.X < r.X

18: L : LLNode

19: R: LLNode

20: end function

21: function Insert(H, s)

22: - inserts a solution s into a high-level tree H

23: H : HLNode

24: s : Solution

25: C ^ new LLNode(s)

26: (G, i) ^ SmallestNonDominatingLayer(H, s) 27: while G = null do

28: Cmin ^ a solution with minimum x from C

29: Cmax ^ a solution with minimum y from C

30: (Tl, T) ^ SPLITX(G.V, Cmin)

31: (Tm, Tr) ^ SplitY(T;, Cmax)

32: G.V ^ Merge(Tl, Merge(C, Tr))

33: if Tm = null then

34: return - no more solutions to push down

35: end if

36: if Tl = null and Tr = null then

37: - the current layer is dominated in whole

38: - just insert pushed solutions as a new layer

39: Insert new HLNode(Tm) after G

40: return

41: end if

42: C ^ Tm

43: G ^ G.N

44: end while

45: Insert new HLNode(C) after last node of H 46: end function

Figure 5: A pseudocode for insertion of a solution into a high-level tree

• If Tm is empty, the remaining layers should remain unchanged. The insertion procedure stops here.

• Otherwise, C ^ Tm , and the insertion procedure continues with the next iteration.

If after the last iterations there are some solutions which were not inserted, a new layer is formed from them and is added as the last layer into the high-level tree.

We omit a long and rigorous proof of correctness of this algorithm for the sake of brevity. The running time of the insertion algorithm sums up from the running time of the lookup algorithm (which is O(log M log(N/ log M))) and

from the total time spent in iterations. Assume that P < M iterations were performed. Without losing generality, assume that the layers of sizes Li... Lp were split in these iterations. Denote the sizes of their pairs after splits to be LL,LM, Lf,..., Lp,LM,Lf. The value LM = 1 corresponds to the initial set C consisting of the solution which is to be inserted. In i-th iteration, the following operations with w(1) complexity were performed:

• finding minimum and maximum of C in O(1+log LMi);

• SplitX in O(1 + log(Lp + LM + Lf));

• SplitY in O(1 + log(LM + Lf));

• inner Merge in O(1 + log(LM-i + Lf));

• outer Merge in O(1 + log(Lp + L^i + Lf)).

In total, the sum of all numbers under logarithms does not exceed 4^p=i Li, and hence is O(N). By Cauchy's inequality, the sum of all running times for all iterations is O(P(1 + log(N/P))). For a fixed N, this function has a maximum when P = ©(N), which both gives us that O(P(1 +log(N/P))) = O(M(1 +log(N/M))) and the worst-case running time of O(N). The layer insertion operations which can happen at the end of the algorithm cost only O(log M) and thus don't change the estimations.

The total running time complexity for the insertion algorithm is:

o(M(l+log^)+logMlogi-^ 2.5 Querying the k-th Solution

The NSGA-II algorithm features a selection operator based on binary tournament which considers every solution exactly twice. It does this by maintaining a permutation of indices 1. . . N, where N is the number of solutions in the population. When a solution is to be taken from the population, it takes an element E of this permutation which was not yet considered and takes a solution which is located at the position of E. When all elements of the permutation are considered, a new (random) permutation is generated instead.

To support this selection operator, we need to support an operation of querying the k-th solution in some predefined order (1 < k < N). Note that querying a random solution can be implemented using this operation (k is chosen randomly uniformly from the range [1; N]). Each query must return a solution, the number of layer it resides at, and the crowding distance of this solution in its layer.

To achieve this aim, we utilize node fields containing subtree sizes, and a field which stores the total size of all keys in the subtree of a high-level tree node. We use the lexicographical order on the solutions: first, the layer number is used to order the solutions, second, the abscissa is used to break ties.

Given k, the number of a solution in the lexicographical order, the first part of the query algorithm finds the layer in which this solution resides. This is done by traversing the high-level tree from the root and tracking the total size of keys in subtrees. The second part uses the same idea to find the solution inside the layer. The third part evaluates the crowding distance and constructs the query result. The algorithm is outlined on Fig. 6. The running time of this algorithm is straightforwardly O(log N).

structure QueryResult

- the result of a query

S : Solution - the solution L : Integer - the layer number C : Real - the crowding distance end structure

function QueryLow(T, k, W, H, L)

- finds k-th solution in a low-level tree T T : LLNode

k:Integer

W : Real - the abscissa difference H : Real - the ordinate difference L : Integer - the layer number if T.L = null then if k < T.L.S then

return QueryLow(T.L, k, W, H, L) else

k — k - T.L.S end if end if

if k = 1 then

R — new QueryResult R.S — T.V, R.L — L if T.P = null or T.N = null then R.C — to

elS6?? \T.N.V.X-T.P.V.X\ |t.,y. :.r.Y-T.P.V.Y\

tt.y W ' H

end if return R else

k—k-1 end if

return QueryLow(T.R, k, W, H, L) end function function Query(T, k)

- returns the k-th solution in a high-level tree T T : HLNode

k:Integer

L : Integer — 0 - the layer number loop

if T.L = null then if k < T.L.W then T — T.L continue else

k — k - T.L.W, L — L + T.L.S end if end if

if k < T.V.S then

LL — leftmost child of T.V RR — rightmost child of T.V W — |RR.V.X - LL.V.X| H — |RR.V.Y - LL.V.Y| return QueryLow(T.V, k, W, H, L) end if

k — k - T.V.S, L — L + 1, T — T.R end loop end function

Figure 6: A pseudocode for querying the k-th solution in lexicographical order

2.6 Deletion of the Worst Solution

In the NSGA-II algorithm, when the number of solutions exceeds the required population size, solutions are deleted from the population starting from the last layer. When one needs to retain some solutions in the last layer, they are deleted in the decreasing order of their crowding distances.

We have not invented an efficient data structure which locates the solution with the smallest crowding distance in a layer. So, to delete a single solution with the worst crowding distance from the last layer, we need to compute crowding distances for every solution in the layer and then delete the one with the smallest crowding distance. This can be done in O(N).

3. EXPERIMENTS

To assess the performance of the proposed implementation, we performed experiments which are largely based on the experiments from [13]. The following benchmark problems are considered: ZDT1-ZDT4 and ZDT6 [17], DTLZ1-DTLZ7 [6], WFG1-WFG9 [8]. The parameters of these problems are the same as in the corresponding papers, except that we considered only two-dimensional versions of these problems.

The NSGA-II algorithm was implemented in as much compatible way to [13] as possible. Specifically, their implementation of the selection operator uses dominance comparison in the first step, unlike classical definition [5], which uses layer number comparison instead. A significant difference is that WFG problems in [13] were implemented using single-precision floating-point numbers, while we used double-precision numbers consistently.

We tested both generational and steady-state variants of the NSGA-II algorithm. Each of these variants was run using the proposed implementation, the ENLU approach [12] and the fast non-dominated sorting from [5]. For each combination, 100 runs were performed, each run had the computational budget of 25 000 evaluations as in [13]. In each run, the random seed for all operations related to evolutionary operators and to individual selection was derived from the number of the run. For other operations which require a random number generator, such as operations with Cartesian trees, a separate random number generator was used.

To assess quality of the results, we evaluated the hypervolume indicator [19] at the end of each run. To assess performance, we evaluated the wall-clock running time and the number of objective comparisons. Notably, each crowding distance evaluation was thought to use four objective comparisons. As the wall-clock time includes not only the time of non-dominated sorting, but the time of evolutionary operators and fitness function evaluations as well, we estimate the latter time by separately running evolutionary operators and fitness functions on randomly generated data and subtracting this time from the wall-clock time. For each mentioned measure we track the median and the interquartile range (IQR).

The source code for experiment reproduction is available both as supplementary materials to the paper and at GitHub1. We ran this code on a server with two Intel® Xeon®E5606 CPUs clocked at 2.13 GHz, having eight cores in total.

1 https: //github.com/mbuzdalov/papers/tree/master/2015-gecco-nsga-ii-ss

For each problem, medians and interquartile ranges for hy-pervolumes were correspondingly equal for all generational algorithms, and the same was true for all steady-state algorithms, so we are not presenting the results for hyper-volumes. All other results of experiments are presented in Table 1.

One can see that the running times and the numbers of comparisons differ significantly. In all cases except one the proposed implementation (called INDS in the tables and further on) performs noticeably faster (in terms of wall-clock time) both with generational and with steady-state versions.

As the running time of Deb's fast non-dominated sorting is not input-sensitive, as follows from the algorithm, the number of comparisons as well as wall-clock times are almost equal for all problems, considering separately generational and steady-state NSGA-II variants. Two other approaches are clearly input-sensitive (for the steady-state case, the biggest time (0.12) is almost three times higher than the smallest time (0.042) in the case of ENLU).

The average operation complexity for INDS is expected to be higher than of the ENLU (due to usage of the tree-based data structure), and ENLU seems to be more complex than the fast non-dominated sorting as well. This made us interested in estimating the overhead imposed by using this or that algorithm, measured as the average wall-clock time per objective comparison. From the data presented in Table 1, we computed this ratio for all three algorithms separately for generational and steady-state variants. It appeared to be almost problem-independent. We give medians and interquartile ranges for these values (notation is the same as in Table 1):

• INDS(gen): median 4.07 - 10-8, IQR 8.40 - 10-10;

• ENLU(gen): median 1.39 - 10-8, IQR 1.90 - 10-10;

• debNDS(gen): median 8.97 - 10-9, IQR 2.31 - 10-11;

• INDS(ss): median 1.36 - 10-8, IQR 1.87 - 10-10;

• ENLU(ss): median 1.08 - 10-8, IQR 2.00 - 10-11;

• debNDS(ss): median 5.88 - 10-9, IQR 1.58 - 10-11.

We account the difference between generational and steady-state variants for different individual removal patterns (bulk removal in the first case and single element removal in the second case). According to the median values of average overhead, INDS is indeed the most expensive (approximately five times more work per single comparison than in fast non-dominated sorting), but it is not significantly worse than ENLU in this aspect. All differences are much smaller for the steady-state variant.

Finally, running times for the same algorithm on the same problem, but for two different NSGA-II variants (generational and steady-state) are very similar in the case of INDS and of ENLU (steady-state times never exceed generational times multiplied by two), but fast non-dominated sorting is slower by more than 10 times in the steady-state case. This removes a big complication of using steady-state multiobjec-tive evolutionary algorithms in practice: while often being superior in terms of quality of results, they are no more slower in terms of running times - at least for biobjective problems.

Table 1: Results of experiments. Notation: "time" stands for the running time (in seconds), "cmp" for the number of comparisons, "med" for the median, "IQR" for the interquartile range. "INDS" is the proposed implementation, "ENLU" is the ENLU-based implementation, "debNDS" is the standard implementation which uses the Deb's fast non-dominated sorting, "gen" corresponds to the generational variant of the NSGA-II algorithm and "ss" to the steady-state variant. Gray cells denote minima of the running time and of the number of comparisons (for generational and steady-state variants separately).

Type INDS (gen) ENLU(gen) debNDS(gen) IND S(ss) ENLU(ss) debNDS (ss)

med IQR med | IQR med | IQR med IQk med | IQR med | IQR

ZDT1

time 5.42 • IO"2 1.13 • 106 1.12 • IO"2 2.68- 104 1.01 • io-1 6.53 • 106 4.47- 10-d 2.98- 104 1.42 • IO"1 1.08 • 10v 2.40 • 10-d 2.17 • 103 7.49 • IO-2 3.70 • 106 2.44- IO-2 7.74 • IO4 9.46 • 10~¿ 7.64- 106 6.11 • io-d 1.21 • 105 1.66 • 10u 2.58- 10s 2.47- 10~¿ 8.06 • IO4

ZDT2

time 5.68- IO-2 1.15 • 106 7.50 • IO-3 2.14- 104 9.16 • 10" ^ 6.10 • 106 9.13 • 10_d 3.95 • 104 1.08 • 10"1 1.07 • 10v 1.78- IO-'1 4.56 • 103 7.09 • IO-2 3.30 • 106 2.53 • IO-3 1.14 • 10s 9.44- 10~¿ 6.90 • 106 4.04 • IO-d 1.61 • 105 1.65 • 10u 2.57- 10s 4.88- IO-'1 1.20 • 10s

ZDT3

time 5.46 • IO-2 1.12 • 106 8.72 • 10-4 2.31 • 104 9.63 • IO-'1 6.32 • 106 5.89 • 10_d 2.63 • 104 9.80 • 10" ^ 1.08 • 10v 1.92 • 10_d 1.76 • 103 7.07- IO-2 3.38- 106 1.80 • IO-3 19.28 • 10~¿ 6.21 • IO4 | 7.14- 106 5.16 • IO-d 8.97 • IO4 1.56 • 10u 2.57- 10s 1.31 • IO-'1 6.77- IO4

ZDT4

time 4.44- IO-2 1.07 • 106 2.89 • IO-3 2.97- 104 5.07 • 3.41 • 106 3.73 • 10_d 3.81 • 104 1.13 • 10"1 1.07 • 10v 3.28- 10_d 8.86 • 103 5.13 • 2.05 3.11 10_d 1.01 • 10f 4.18- 10" 2 3.68- 106 4.48 • IO-3 1.60 • 10s 1.78- 10u 2.56 • 10s 5.241.20 • 105

ZDT6

time 5.61 • IO-2 1.12 • 106 2.12 • IO-3 1.11 • 104 6.80 • 10" ^ 4.52 • 106 6.43 • 10-4 3.89 • 104 1.20 • 10"1 1.07 • 10v 2.51 • 10_d 4.40 • 103 6.65 • IO-2 2.64- 106 2.35 • IO-3 9.39 • IO4 7.19 • 10" ^ 5.29 • 106 2.91 • 10_d 1.13 • 105 1.71 • 10u 2.56 • 10s 3.26 • 9.17- IO4

DTLZ1

time 4.67- IO-2 1.02 • 106 5.04 • IO-3 2.27- 104 5.66 • 3.73 • 106 3.19 • 10_d 5.80 • 104 1.13 • 10"1 1.07 • 10v 2.86 • 10_d 7.87 • 103 4.88- IO-2 2.21 • 106 3.08- IO-3 2.02 • 105 5.21 • 4.17- 106 5.21 • 10_d 2.42 • 105 1.69 • 10u 2.56 • 10s 2.83 • 2.19 • 105

DTLZ2

time 5.41 • 10~2 1.09 • 106 4.08 • 10~3 2.81 • 104 9.51 • 10~¿ 6.24 • 106 4.11 • 10-d 2.29 • 104 9.71 • 10~¿ 1.08 • 10v 3.78- 10-d 1.24 • 103 6.75 • IO-2 4.19 • 106 3.92 • IO-3 5.00 • IO4 9.35 • IO-^ 7.91 • 106 5.37 • IO" d 1.77 • 105 1.53 • 10u 2.58- 10s 2.15 • 10~¿ 4.97- IO4

DTLZ3

time 4.30 • IO"2 1.08 • 106 1.83 • 10~3 3.78- 104 5.37 • 10" ^ 3.81 • 106 4.27- 10-d 1.11 • 105 1.13 • IO"1 1.06 • 10v 9.16 • IO-4 8.25 • 103 4.18- IO-2 1.31 • 106 2.07- 10~3 1.06 • 105 4.78- IO-^ 3.19 • 106 4.52 • IO" d 1.02 • 105 1.94- 10u 2.55 • 10s 5.69 • 10~¿ 1.07• 105

DTLZ4

time 5.60 • IO-2 1.08 • 106 2.81 • IO-3 3.32 • 104 9.47 • 10" ^ 6.21 • 106 3.34- 10_d 3.57- 104 1.02 • 10"1 1.08 • 10v 9.92 • 10_d 4.89 • 103 6.87- IO-2 3.95 • 106 6.93 • IO-3 8.01 • IO4 9.30 • IO-^ 7.67- 106 4.99 • IO-d 2.31 • 105 1.58- 10u 2.58- 10s 5.08- IO-'1 8.87 • IO4

DTLZ5

time 5.36 • IO-2 1.09 • 106 2.49 • 10-4 3.29 • 104 9.30 • 10" ^ 6.24 • 106 1.24- 10_d 2.29 • 104 9.66 • IO-'1 1.08 • 10v 3.03 • 10_d 1.24 • 103 6.76 • IO-2 4.19 • 106 4.36 • IO-3 19.40 • 10~A 5.26 • IO4 | 7.91 • 106 4.35 • IO-d 1.77 • 105 1.54- 10u 2.58- 10s 2.02 • IO-'1 4.97- IO4

DTLZ6

time 4.96 • IO-2 1.16 • 106 6.55 • IO-3 2.38- 104 7.09 • 10" ^ 5.17 • 106 1.94- 10_d 4.81 • 104 1.03 • 10"1 1.07 • 10v 1.79 • 10_d 4.27 • 103 5.88- IO-2 2.90 • 106 2.52 • IO-3 2.13 • 10s 7.18- IO-^ 5.95 • 106 5.43 • 10_d 2.46 • 10s 1.70 • 10u 2.57- 10s 1.96 • IO-'1 2.13 • 105

DTLZ7

time 5.64- IO-2 1.13 • 106 4.39 • IO-3 2.07- 104 9.43 • 6.12 • 106 3.82 • 10_d 2.86 • 104 1.02 • 10"1 1.08 • 10v 1.80 • 10-d 1.85 • 103 7.06 • IO-2 3.55 • 106 2.16 • IO-3 6.15 • IO4 9.077.21 • 106 2.93 • 10_d 1.32 • 105 1.59 • 10u 2.57- 10s 1.65 • 6.50 • IO4

WFG1

time 4.90 • IO-2 1.07 • 106 1.48 • IO-3 1.99 • 104 7.91 • 5.48 • 106 1.87- 10_d 1.62 • 105 1.00 • 10"1 1.07 • 10v 1.12 • 10-d 1.14 • 104 6.11 • IO-2 2.94- 106 4.02 • IO-3 4.30 • 105 8.15 • IO-^ 6.27- 106 1.15 • IO-^ 4.71 • 105 1.57- 10u 2.57- 10s 2.80 • 5.18- 105

WFG2

time 5.94- 10~2 1.07 • 106 2.67 • IO-3 3.65 • 104 1.04 • IO"1 6.97 • 106 5.11 • 10~3 6.68- 104 1.02 • IO"1 1.08 • 10v 1.87- 10~3 2.14 • 103 9.73 • IO-2 7.05 • 106 4.72 • IO-3 1.25 • 105 1.23 • IO"1 1.13 • 10v 8.75 • IO" 3 1.64 • 105 1.54- 10u 2.61 • 10s 1.98- IO-2 8.74 • IO4

WFG3

time 5.63 • IO"2 1.07 • 106 3.37 • 10~3 2.69 • 104 9.40 • 10" ^ 6.16 • 106 5.75 • 10~4 2.88- 104 9.83 • 10~¿ 1.08 • 10v 2.29 • 10~d 1.42 • 103 7.59 • IO-2 4.98- 106 5.90 • IO-3 7.89 • IO4 9.56 • 10~¿ 8.69 • 106 6.60 • io-d 1.94 • 105 1.53 • 10u 2.59 • 10s 1.72 • 10~¿ 7.67- IO4

WFG4

time 5.34- IO"2 1.07 • 106 1.46 • 10~3 2.11 • 104 9.22 • 10~¿ 6.07 • 106 2.35 • 10~d 5.90 • 104 9.90 • 10~¿ 1.08 • 10v 1.19 • 10_d 2.56 • 103 7.51 • IO-2 4.76 • 106 4.20 • IO-3 1.92 • 105 9.78- 10~ ¿ 8.43 • 106 3.68 • 10-d 2.49 • 105 1.53 • 10u 2.59 • 10s 2.04- 10~¿ 1.85 • 105

WFG5

time 5.62 • IO"2 1.08 • 106 4.48 • 10-3 4.42 • 104 1.03 • io-1 7.10 • 106 2.08- 10_d 1.80 • 10s 9.73 • 10~¿ 1.08 • 10v 5.45 • 10_d 7.08 • 103 9.01 • IO-2 6.71 • 106 3.21 • 10~3 4.86 • 105 1.20 • IO"1 1.11 • 10v 3.73 • 10-d 7.35 • 10s 1.54- 10u 2.61 • 10s 1.60 • 10~¿ 4.73 • 10s

WFG6

time 5.46 • IO"2 1.07 • 106 3.50 • 10~3 1.97- 104 9.19 • IO-^ 5.96 • 106 6.80 • 10~d 7.38- 104 1.01 • IO"1 1.08 • 10v 3.19 • 10~d 4.13 • 103 6.92 • IO-2 4.48- 106 2.74- 10~3 3.28 • 105 9.29 • 10~¿ 8.10 • 106 7.66 • IO" d 4.59 • 10s 1.53 • 10u 2.58- 10s 1.88- 10~¿ 3.29 • 10s

WFG7

time 5.58- IO-2 1.06 • 106 1.64 • IO-3 2.98- 104 9.80 • 10_ ^ 6.38 • 106 2.45 • 10_d 2.43 • 104 9.96 • 10"^ 1.08 • 10v 4.25 • 10-d 1.11 • 103 7.81 • IO-2 5.53 • 106 1.67- IO-3 7.61 • IO4 1.02 • IO-1 9.47- 106 2.84 • IO-d 1.35 • 10s 1.52 • 10u 2.59 • 10s 1.63 • 7.71 • IO4

WFG8

time 5.50 • IO-2 1.05 • 106 4.76 • 10"3 4.14- 104 8.04 • IO-^ 4.59 • 106 1.46 • 10_d 1.25 • 10s 1.08 • IO-1 1.07 • 10v 1.40 • 10-d 8.61 • 103 4.79 • IO-2 1.66 • 106 3.04- IO-3 1.32 • 105 6.284.12 • 106 6.48 • IO-d 1.61 • 10s 1.52 • 10u 2.56 • 10s 2.641.23 • 10s

WFG9

!mp 5.67- 10~2 1.06 • 106 2.99 • 10~3 2.69 • 104 9.68 • IO-2 6.23 • 106 3.83 • 10~3 1.05 • 10s 1.01 • io-1 1.08 • 10v 5.63 • 10~3 4.66 • 103 7.46 • IO-2 4.91 • 106 3.45 • IO-3 4.58 • 105 1.04- IO"1 8.56 • 106 3.93 • IO-3 6.21 • 10s 1.53 • 10u 2.59 • 10s 1.56 • IO-2 4.31 • 10s

4. CONCLUSION

We proposed a new approach to implementation of steady-state multiobjective evolutionary algorithms for two dimensions. This approach is based on a data structure which is able to perform fast incremental non-dominated sorting and track other valuable data such as crowding distance. Experiments with the NSGA-II algorithm showed that the new approach offers running times similar to or better than those of generational versions while retaining the quality of steady-state versions.

This work was financially supported by the Government of Russian Federation, Grant 074-U01.

5. REFERENCES

[1] H. A. Abbass, R. Sarker, and C. Newton. PDE: A Pareto Frontier Differential Evolution Approach for Multiobjective Optimization Problems. In Proceedings of the Congress on Evolutionary Computation, pages 971-978. IEEE Press, 2001.

[2] M. Buzdalov and A. Shalyto. A provably asymptotically fast version of the generalized Jensen algorithm for non-dominated sorting. In International Conference on Parallel Problem Solving from Nature, number 8672 in Lecture Notes in Computer Science, pages 528-537. 2014.

[3] D. W. Corne, N. R. Jerram, J. D. Knowles, and M. J. Oates. PESA-II: Region-based Selection in Evolutionary Multiobjective Optimization. In Proceedings of Genetic and Evolutionary Computation Conference, pages 283-290. Morgan Kaufmann Publishers, 2001.

[4] D. W. Corne, J. D. Knowles, and M. J. Oates. The Pareto Envelope-based Selection Algorithm for Multiobjective Optimization. In Parallel Problem Solving from Nature Parallel Problem Solving from Nature VI, number 1917 in Lecture Notes in Computer Science, pages 839-848. Springer, 2000.

[5] K. Deb, A. Pratap, S. Agarwal, and T. Meyarivan. A Fast Elitist Multi-Objective Genetic Algorithm: NSGA-II. Transactions on Evolutionary Computation, 6:182-197, 2000.

[6] K. Deb, L. Thiele, M. Laumanns, and E. Zitzler. Scalable Test Problems for Evolutionary Multiobjective Optimization, pages 105-145. Springer, 2005.

[7] F.-A. Fortin, S. Grenier, and M. Parizeau. Generalizing the Improved Run-time Complexity Algorithm for Non-dominated Sorting. In Proceeding of the Fifteenth Annual Conference on Genetic and

Evolutionary Computation Conference, GECCO '13, pages 615-622. ACM, 2013.

[8] S. Huband, P. Hingston, L. Barone, and R. L. While. A review of multiobjective test problems and a scalable test problem toolkit. IEEE Transactions on Evolutionary Computation, 10(5):477-506, 2006.

[9] M. T. Jensen. Reducing the Run-time Complexity of Multiobjective EAs: The NSGA-II and Other Algorithms. Transactions on Evolutionary Computation, 7(5):503-515, 2003.

[10] J. D. Knowles and D. W. Corne. Approximating the Nondominated Front Using the Pareto Archived Evolution Strategy. Evolutionary Computation, 8(2):149-172, 2000.

[11] H. T. Kung, F. Luccio, and F. P. Preparata. On finding the maxima of a set of vectors. Journal of ACM, 22(4):469-476, 1975.

[12] K. Li, K. Deb, Q. Zhang, and S. Kwong. Efficient non-domination level update approach for steady-state evolutionary multiobjective optimization. Technical report, 2014.

[13] A. J. Nebro and J. J. Durillo. On the effect of applying a steady-state selection scheme in the multi-objective genetic algorithm NSGA-II. In Nature-Inspired Algorithms for Optimisation, number 193 in Studies in Computational Intelligence, pages 435-456. Springer Berlin Heidelberg, 2009.