Реализация на C# алгоритма Флойда-Уоршелла нахождения всех кратчайших путей между всеми вершинами графа

2
826
views

Вступление

Задача поиска кратчайшего пути на графе встречается повсеместно. Мы сталкиваемся с ней, когда идем на работу в ISsoft, когда открываем браузер или, например, когда раскладываем предметы на своем столе поудобнее.

Звучит капельку чересчур, правда?

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

В приведенных примерах «короткость» пути определяется или преодолеваемым расстоянием, или затрачиваемым интервалом времени. В примере со столом «короткость» пути может определяться количеством (или сложностью) необходимых действий: открыть ящик стола, открыть папку, достать лист бумаги – против взять лист из стопки на столе.

Если внимательно посмотреть, то каждый из приведенных примеров описывает поиск кратчайшего пути из одной исходной точки. Такой класс задач называется SSSP (Single Source Shortest Path) и базовым алгоритмом решения таких задач, является алгоритм Дейкстры, имеющий квадратичную сложность. Именно такие задачи мы решаем чаще всего. 

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

  1. Работа — Дом
  2. Дом — Работа (пути могут отличаться, например из-за наличия улиц с односторонним движением)
  3. Работа — Театр
  4. Театр — Работа
  5. Дом — Театр
  6. Театр — Дом

С увеличением интересующих нас мест количество пар маршрутов, будет стремительно возрастать, согласно формулам комбинаторики, если точнее, то формуле размещения:

A(k, n) = n! / (n - m)!
// где n - количество элементов, 
//     k - количество элементов в размещение (в нашем случае k = 2)

Не трудно посчитать, что между 4-я местами существует 12 маршрутов, а между 10-ю их уже 90. И все это не учитывая факта, что между интересующими нас путями может не быть прямого пути, а это привносит в нашу схему промежуточные места (напр. остановки, перекрестки). Класс задач, описывающих поиск всех кратчайших путей между всеми вершинами графа (каждое место представляет собой вершину графа) называется APSP (All Pairs Shortest Paths) и базовым алгоритмом решения таких задач является алгоритм Флойда-Уоршелла, которые мы и рассмотрим.

Алгоритм Флойда-Уоршелла

Алгоритм Флойда-Уоршелла нахождения всех кратчайших путей между всеми вершинами графа был представлен Робертом Флойдом в работе [1] (см. секцию «Ссылки» в самом конце). В том же году, Питер Ингерман в работе [2] представил современную, программную формулировку алгоритма в виде трех вложенных циклов for, приблизительный псевдокод которой приведен ниже:

algorithm FloydWarshall(W) do
  for k = 0 to N - 1 do
    for i = 0 to N - 1 do 
      for j = 0 to N - 1 do 
        W[i,j] = min(W[i,j], W[i,k] + W[k,j])
	end for
    end for
  end for
end algorithm

// где W     - матрица весов размером N x N, 
//     min() - функция возвращающая минимальный из переданных аргументов.

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

На рисунке ниже изображен взвешенный, ориентированный граф из 5-и вершин. Граф изображен как графически (слева), так и в виде матрицы весов (справа). Матрица весов представляет собой вариант матрицы смежности, где каждая ячейка матрицы содержит «вес» — расстояние из вершины i (номер строки) в вершину j (номер столбца). Важно отметить, что матрица весов не содержит «путь» (список вершин, через которые пролегает путь между вершинами) — только расстояние между ними.

На рисунке все ячейки матрицы W содержат значения, равные весам ребер между смежными вершинами графа. Поэтому, если изучить пути из вершины 0 (строка 0), то окажется, что их … нет. Однако, на графическом представлении графа отчетливо видно наличие таких путей, правда, проходящих через другие вершины графа — из вершины 0, через вершину 1, можно попасть в вершины 2 и 3. Причина этого отличия в том, что в изначальном состоянии, матрица весов содержит исключительно информацию о расстояниях между смежными вершинами графа.

Однако представленной в матрице информации на самом деле достаточно, чтобы обнаружить все остальные пути. Обратите внимание на ячейку W[0,1]. Значение в этой ячейке говорит о существовании пути 0→1:1 (из вершины 0 в вершину 1, с весом 1). Наличие этого пути означает, что мы можем просканировать все ячейки строки 1 матрицы весов (в которой содержатся все пути из вершины 1) и добавить эту информацию в строку 0, предварительно увеличив их на вес записанный в ячейке W[0,1] (на расстояние из вершины 0 в вершину 1):

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

На рисунке отчетливо видно, что известный нам путь 0→3 имеет вес 4 (помним, что матрица весов не содержит «самого пути» т.е. мы не знаем, через какие вершины проходит этот путь), в то время, как новый путь 0→2→3 имеет вес 3. В этот момент мы совершаем выбор и записываем в ячейку W[0,3] меньшее значение — 3.

Кажется, мы только что нашли кратчайший путь из вершины 0 в вершину 3, ну или хотя бы один из них 🙂

Именно так и работает алгоритм Флойда-Уоршелла. Посмотрите на псевдокод алгоритма еще раз:

algorithm FloydWarshall(W) do
  for k = 0 to N - 1 do
    for i = 0 to N - 1 do 
      for j = 0 to N - 1 do 
	  W[i,j] = min(W[i,j], W[i,k] + W[k,j])
	end for
    end for
  end for
end algorithm

Внешний цикл for по k обходит все вершины графа. k представляет собой текущую вершину, через которую мы выполняем поиск. Внутренний цикл for по i обходит все вершины графа, из которых мы выполняем поиск путей, проходящих через вершину k. И наконец внутренний цикл for по j обходит все вершины графа, в которые мы выполняем поиск из вершины i через вершину k. И наконец, в теле цикла, мы сравниваем вес известного нам пути i→j записанного в матрице весов с весом пути i→k→j, записывая в ячейку минимальное из значений. 

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

Базовая реализация алгоритма Флойда-Уоршелла

Внимание! Программный код приведенных ниже реализаций алгоритма и используемые в экспериментах графы можно найти в репозитории на GitHub.

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

  1. Все приведенные реализации алгоритма работают с матрицей весов, представленной в виде массива.
  2. Все приведенные реализации алгоритма используют целочисленную арифметику. Отсутствие пути между вершинами представляется в матрице константой NO_EDGE = (int.MaxValue / 2) - 1.

Теперь по прядку. Первый пункт.

Когда мы говорим о матрицах, мы оперируем такими понятиями как «строка» и «столбец». Из-за этого нам естественно думать о матрице как о «квадрате» или «прямоугольнике», как показано в подпункте а) на рисунке ниже. 

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

i = row * row_size + col;
// где row      - индекс строки, в которой расположена ячейка,
//     col      - индекс столбца, в котором расположена ячейка, 
//     row_size - количество ячеек матрицы в одной строке.

Представление матрицы в виде массива является необходимым предусловием для эффективного выполнения алгоритма Флойда-Уоршелла. И обусловлена она тем, что такое представление повышает локализацию доступа к данным. Тема эффективного доступа к памяти невероятно объемная и является объектом многих научных исследований в области разработки алгоритмов. Её просто невозможно описать в нескольких абзацах, поэтому в этот раз я попрошу вас поверить мне на слово (а перед этим, разумеется, изучить тему самостоятельно: не стоит доверять всяким незнакомцам в интернете).

Насчет второго пункта.

Если внимательно посмотреть, на псевдокод алгоритма, можно заметить, что там нигде нет проверки, что путь между двумя вершинами существует. Вместо этого в псевдокоде используется функция min(), возвращающая минимальное значение. Причина в том, что если между двумя вершинами в графе нет пути, то ячейка матрицы весов должна содержать в себе значение бесконечности (а любое число, не являющееся бесконечностью, меньше бесконечности — если это конечно не JavaScript, там возможно все). При использовании целых чисел может появиться соблазн использовать максимальное допустимое значение — int.MaxValue. Однако в этом случае, мы столкнемся с переполнением, когда веса путей i→k и k→j будут равны int.MaxValue (при использовании константы (int.MaxValue / 2) - 1 переполнения не будет). 

«Окей, но ведь можно просто предварительно проверять, существует ли пути i→k и k→j, например сравнивая их с 0 перед выполнением».

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

Уф, с теорией закончили — теперь к реализации:

public void FloydWarshall_00(int[] matrix, int sz) 
{
  for (var k = 0; k < sz; ++k)
  { 
    for (var i = 0; i < sz; ++i)
    {
      for (var j = 0; j < sz; ++j)
      {
        var distance = matrix[i * sz + k] + matrix[k * sz + j];
        if (matrix[i * sz + j] > distance)
        {
          matrix[i * sz + j] = distance;
        }
      }
    }
  }
}

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

Пришло время запустить код на выполнение! Посмотрим, насколько быстро он работает! 

То, что он работает, можно убедиться запустив тесты в репозитории.

В качестве платформы для запуска я буду использовать Benchmark.NET (версия 0.12.1). Экспериментальные графы представляют собой ориентированные, ациклические графы (размером 300, 600, 1200, 2400 и 4800 вершин). Количество ребер в графах приближается к 80% от возможного (максимальное количество ребер в ориентированном, ациклическом графе вычисляется по формуле: (v * (v - 1)) / 2, где v — это количество вершин в графе).

Поехали!

Результаты выполнения на моем компьютере (ОС Windows 10.0.19042, процессор Intel Core i7-7700 CPU 3.60GHz (Kaby Lake) / 8 логических процессоров / 4 физических ядра) приведены в таблице ниже:

MethodSizeMeanErrorStdDev
FloydWarshall_00 300 27.525 ms0.1937 ms0.1617 ms
FloydWarshall_00 600 217.897 ms1.6415 ms1.5355 ms
FloydWarshall_00 1200 1,763.335 ms7.4561 ms6.2262 ms
FloydWarshall_00 2400 14,533.335 ms63.3518 ms52.9016 ms
FloydWarshall_00 4800 119,768.219 ms181.5514 ms160.9406 ms


Из полученных результатов примечательно, что время выполнения невероятно возрастает с увеличением размера графа — так поиск всех кратчайших путей между всеми вершинами в графа из 300 вершин занимает 27 миллисекунд, в графе из 2400 вершин — 14.5 секунд, а в графе из 4800 — целых 119 секунд! Почти 2 минуты!

Смотря на реализацию алгоритма, трудно представить, что мы можем что-то изменить и сократить время выполнения алгоритма. Там ведь всего ТРИ цикла! 

Однако, как это очень часто бывает, возможности скрываются в деталях 🙂

Знай свои данные — не плотные графы

В базовой формулировке алгоритм Флойда-Уоршелла является наиболее эффективным при использовании плотных или полных графов т. к. по своей механике во время выполнения алгоритм обходит каждую вершину графа через каждую вершину. Однако, используемые нами графы, хоть и являются «достаточно плотными» (помним про 80%), в первую очередь являются ациклическими. Простыми словами это означает, что если из вершины 1 есть путь в вершину 2, то из вершины 2 не существует пути в вершину 1. А это в свою очередь означает, наличие значительного количества вершин между которыми нет и не может быть пути. Мы можем это использовать и не проверять наличие путей из вершины i в вершины j, если отсутствует путь i→k:

public void FloydWarshall_01(int[] matrix, int sz)
{
  for (var k = 0; k < sz; ++k)
  {
    for (var i = 0; i < sz; ++i)
    {
      if (matrix[i * sz + k] == NO_EDGE)
      {
        continue;
      }
      for (var j = 0; j < sz; ++j)
      {
        var distance = matrix[i * sz + k] + matrix[k * sz + j];
        if (matrix[i * sz + j] > distance)
        {
          matrix[i * sz + j] = distance;
        }
      }
    }
  }
}

Результаты выполнения предыдущей (00) обновленной (01) реализации приведены в таблице ниже:

MethodSizeMeanErrorStdDevRatio
FloydWarshall_00 30027.525 ms0.1937 ms0.1617 ms1.00
FloydWarshall_01 30012.399 ms0.0943 ms0.0882 ms0.45
FloydWarshall_00 600217.897 ms1.6415 ms1.5355 ms1.00
FloydWarshall_01 60099.122 ms0.8230 ms0.7698 ms0.45
FloydWarshall_00 12001,763.335 ms7.4561 ms6.2262 ms1.00
FloydWarshall_01 1200766.675 ms6.1147 ms5.7197 ms0.43
FloydWarshall_00 240014,533.335 ms63.3518 ms52.9016 ms1.00
FloydWarshall_01 24006,507.878 ms28.2317 ms26.4079 ms0.45
FloydWarshall_00 4800119,768.219 ms181.5514 ms160.9406 ms1.00
FloydWarshall_01 480055,590.374 ms414.6051 ms387.8218 ms0.46


Выглядит впечатляюще! 

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

Но можно ли сделать что-нибудь еще? 🙂

Знай свое оборудование — параллелизм

Мой компьютер оборудован процессором Intel Core i7-7700 CPU 3.60GHz (Kaby Lake) с 8-ю логическими процессорами (HW) или 4-я физическими ядрами с технологий Hyper Threading. Наличие нескольких ядер означает, что мы имеем в своем распоряжении «свободные рабочие руки», которые можно (и нужно) занять полезной работой. Однако, прежде чем создавать потоки, нам необходимо понять — какую именно часть работы алгоритма мы можем разделить.

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

Начнем с цикла по k. На каждой итерации цикла алгоритм рассчитывает пути из каждой вершины в каждую, через вершину k, записывая обнаруженные пути в матрицу весов. По отдельности каждая из итераций выглядит независимой. Например, не имеет значения в каком порядке выполняются итерации k-го цикла: 0, 1, 2, 3 или 3, 1, 2, 0. Однако, основным условием является то, что все эти итерации должны выполняться последовательно (одна за другой). Только тогда выполнится основное условие алгоритма: обход графа из каждой вершины в каждую вершину через каждую вершину.

Теперь рассмотрим цикл по i. На каждой итерации мы ищем пути из вершины i в вершину k и затем из вершины k в вершину j. Если обратиться к матрице весов, то можно заметить, что на каждой итерации цикла мы считываем значение с k-й строки матрицы и обновляем значения только в i-й строке матрицы. Это означает, что итерации цикла по i не только могут быть выполнены в любом порядке (как и итерации цикла по k), но и не используют результат работы друг друга, а значит, могут выполняться параллельно!

Кажется, мы что-то нашли — реализуем это!

Цикл по j, так же как и цикл по i может быть распараллелен. Однако, распараллеливание вложенного цикла является крайне неэффективным. Чтобы не брать мои слова на веру, вы можете проверить это сами используя программный код из репозитория.

public void FloydWarshall_02(int[] matrix, int sz)
{
  for (var k = 0; k < sz; ++k)
  {
    Parallel.For(0, sz, i =>
    {
      if (matrix[i * sz + k] == NO_EDGE)
      {
        return;
      }
      for (var j = 0; j < sz; ++j)
      {
        var distance = matrix[i * sz + k] + matrix[k * sz + j];
        if (matrix[i * sz + j] > distance)
        {
          matrix[i * sz + j] = distance;
        }
      }
    });
  }
}

В нашей реализации распараллеливание выполняется при помощи класса Parallel. Результаты выполнения алгоритмов 01 и 02 приведены в таблице ниже (колонка Ratio показывает значение относительно алгоритма 00):

MethodSizeMeanErrorStdDevRatio
FloydWarshall_01 30012.399 ms0.0943 ms0.0882 ms0.45
FloydWarshall_02 3006.252 ms0.0211 ms0.0187 ms0.23
FloydWarshall_01 60099.122 ms0.8230 ms0.7698 ms0.45
FloydWarshall_02 60035.825 ms0.1003 ms0.0837 ms0.16
FloydWarshall_01 1200766.675 ms6.1147 ms5.7197 ms0.43
FloydWarshall_02 1200248.801 ms0.6040 ms0.5043 ms0.14
FloydWarshall_01 24006,507.878 ms28.2317 ms26.4079 ms0.45
FloydWarshall_02 24002,076.403 ms20.8320 ms19.4863 ms0.14
FloydWarshall_01 480055,590.374 ms414.6051 ms387.8218 ms0.46
FloydWarshall_02480015,614.506 ms115.6996 ms102.5647 ms0.13


Из полученных результатов отчетливо видно — распараллеливание сократило время выполнения алгоритма от 2-х до 4-х раз! Однако важно помнить, что выполнение параллельной реализации потребляет практически все доступные на компьютере ресурсы процессора, что может привести к их серьезному дефициту у других приложений в системе.

Возможно вы догадались — это еще не все 🙂

Знай свою платформу — векторизация

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

Звучит сложнее, чем оно есть на самом деле, поэтому давайте разберемся на примере:

var a = new [] { 5, 7, 8, 1 };
var b = new [] { 4, 2, 2, 6 };
var c = new [] { 0, 0, 0, 0 };

for (var i = 0; i < 4; ++i)
  c[i] = a[i] + b[i];

В сверх упрощённом виде итерация 0 цикла for выглядит следующим образом:

Первым шагом – значения элементов массивов a и b из памяти загружаются в регистры процессора (обратите внимание, что в один регистр помещается только один элемент массива). Как только значения оказались в регистрах, над ними выполняется операция сложения и ее результат записывается обратно в память в соответствующий элемент массива c. Подобным образом выполняются остальные итерации цикла, т. е. этот процесс повторяется 4 раза.

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

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

В .NET использование векторных инструкций и векторных регистров доступно через специальные типы — Vector и Vector<T> (а также, через типы из пространства имен System.Runtime.Intrinsics). Реализация алгоритма с использованием типов Vector приведена ниже:

public void FloydWarshall_03(int[] matrix, int sz)
{
  for (var k = 0; k < sz; ++k)
  {
    for (var i = 0; i < sz; ++i)
    {
      if (matrix[i * sz + k] == NO_EDGE)
      {
        continue;
      }

      var ik_vec = new Vector<int>(matrix[i * sz + k]);

      var j = 0;
      for (; j < sz - Vector<int>.Count; j += Vector<int>.Count)
      {
        var ij_vec  = new Vector<int>(matrix, i * sz + j);
        var ikj_vec = new Vector<int>(matrix, k * sz + j) + ik_vec;

        var lt_vec = Vector.LessThan(ij_vec, ikj_vec);
        if (lt_vec == new Vector<int>(-1))
        {
          continue;
        }

        var r_vec = Vector.ConditionalSelect(lt_vec, ij_vec, ikj_vec);
        r_vec.CopyTo(matrix, i * sz + j);
      }

      for (; j < sz; ++j)
      {
        var distance = matrix[i * sz + k] + matrix[k * sz + j];
        if (matrix[i * sz + j] > distance)
        {
          matrix[i * sz + j] = distance;
        }
      }
    }
  }
}

Первым шагом var ik_vec = new Vector<int>(matrix[i * sz + k]) в векторизованной части программного кода мы создаем новый вектор, в котором каждый элемент равен весу пути i→k. В результате, если в векторе можно разместить 4-е значения типа int, а вес пути i→k равен 5, то созданный вектор будет иметь значение [5,5,5,5]. Затем, на каждой итерации мы одновременно рассматриваем пути из вершины i в вершины j,j + 1,...,j + Vector<int>.Count.

Использование Vector<int>.Count возвращает возможное количество элементов типа int размещаемых векторном регистре. В зависимости от модели процессора векторные регистры могут иметь разный размер.

Важное отличие векторной реализации в том, что мы работаем над несколькими элементами одновременно, и поэтому результат операции может быть не однозначен. Например, операция сравнения Vector.LessThan(ij_vec,ikj_vec) возвращает новый вектор, каждый элемент которого содержит результат сравнения соответствующих элементов из векторов ij_vec и ikj_vec (-1, если элемент из ij_vec меньше элемента из ikj_vec и 0, если наоборот). В свою очередь, операция Vector.ConditionalSelect(lt_vec,ij_vec,ikj_vec) возвращает новый вектор, который содержит минимальные элементы из обоих векторов, т.е. если элемент вектора lt_vec установлен в -1, то берется элемент из ij_vec, иначе выбирается соответствующий элемент из ikj_vec. И наконец, операция r_vec.CopyTo(matrix,i * sz + j) записывает значение вектора обратно в массив.

Разбиение реализации на 2 цикла (векторизованный и не векторизованный) необходимо из-за случаев, когда длина массива не кратна значению возвращаемому Vector<int>.Count. В этом случае мы не можем обработать все элементы массива, используя векторные инструкции, и нам нужно рассчитать оставшиеся пути, использовав обычные инструкции.

Результатом наших трудов стали следующие результаты (сравнение реализации 01, 02 и 03):

MethodSizeMeanErrorStdDevRatio
FloydWarshall_01 30012.399 ms0.0943 ms0.0882 ms0.45
FloydWarshall_02 3006.252 ms0.0211 ms0.0187 ms0.23
FloydWarshall_033003.056 ms0.0301 ms0.0281 ms0.11
FloydWarshall_01 60099.122 ms0.8230 ms0.7698 ms0.45
FloydWarshall_02 60035.825 ms0.1003 ms0.0837 ms0.16
FloydWarshall_0360024.378 ms0.4320 ms0.4041 ms0.11
FloydWarshall_01 1200766.675 ms6.1147 ms5.7197 ms0.43
FloydWarshall_02 1200248.801 ms0.6040 ms0.5043 ms0.14
FloydWarshall_031200185.628 ms2.1240 ms1.9868 ms0.11
FloydWarshall_01 24006,507.878 ms28.2317 ms26.4079 ms0.45
FloydWarshall_02 24002,076.403 ms20.8320 ms19.4863 ms0.14
FloydWarshall_0324002,568.676 ms31.7359 ms29.6858 ms0.18
FloydWarshall_01 480055,590.374 ms414.6051 ms387.8218 ms0.46
FloydWarshall_02480015,614.506 ms115.6996 ms102.5647 ms0.13
FloydWarshall_03480018,257.991 ms84.5978 ms79.1329 ms0.15


Из полученных результатов видно, что использование векторных инструкций в значительной степени сократило время выполнение алгоритма: от 3-х до 4-х раз по сравнению предыдущей не использующей параллелизм версией (01). Интересно, что векторизованная версия демонстрирует преимущество над распараллеленной версией на небольших размерах графов (до 2400 вершин).

И напоследок — совместим векторизацию и параллелизм!

Если вас заинтересовал вопрос векторизации, то я рекомендую ознакомиться с циклом статей от Dan Shechter, в которых он занимается векторизацией метода Array.Sort (результаты работы представленные в его статьях впоследствии нашли свое место в реализации сборщика мусора в .NET 5).

Знаю свою платформу и свое оборудование — векторизация и параллелизм!

Последняя реализация совмещает в себе использование параллелизма и векторизации и … она не отличается оригинальностью 🙂 По сути, мы просто вставляем тело цикла по i из реализации 03 в тело метода Parallel.For из реализации 02:

public void FloydWarshall_04(int[] matrix, int sz)
{
  for (var k = 0; k < sz; ++k)
  {
    Parallel.For(0, sz, i =>
    {
      if (matrix[i * sz + k] == NO_EDGE)
      {
        return;
      }

      var ik_vec = new Vector<int>(matrix[i * sz + k]);

      var j = 0;
      for (; j < sz - Vector<int>.Count; j += Vector<int>.Count)
      {
        var ij_vec  = new Vector<int>(matrix, i * sz + j);
        var ikj_vec = new Vector<int>(matrix, k * sz + j) + ik_vec;

        var lt_vec = Vector.LessThan(ij_vec, ikj_vec);
        if (lt_vec == new Vector<int>(-1))
        {
          continue;
        }

        var r_vec = Vector.ConditionalSelect(lt_vec, ij_vec, ikj_vec);
        r_vec.CopyTo(matrix, i * sz + j);
      }

      for (; j < sz; ++j)
      {
        var distance = matrix[i * sz + k] + matrix[k * sz + j];
        if (matrix[i * sz + j] > distance)
        {
          matrix[i * sz + j] = distance;
        }
      }
    });
  }
}

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

MethodSizeMeanErrorStdDevRatio
FloydWarshall_0030027.525 ms0.1937 ms0.1617 ms1.00
FloydWarshall_01 30012.399 ms0.0943 ms0.0882 ms0.45
FloydWarshall_02 3006.252 ms0.0211 ms0.0187 ms0.23
FloydWarshall_033003.056 ms0.0301 ms0.0281 ms0.11
FloydWarshall_00600217.897 ms1.6415 ms1.5355 ms1.00
FloydWarshall_01 60099.122 ms0.8230 ms0.7698 ms0.45
FloydWarshall_02 60035.825 ms0.1003 ms0.0837 ms0.16
FloydWarshall_0360024.378 ms0.4320 ms0.4041 ms0.11
FloydWarshall_0012001,763.335 ms7.4561 ms6.2262 ms1.00
FloydWarshall_01 1200766.675 ms6.1147 ms5.7197 ms0.43
FloydWarshall_02 1200248.801 ms0.6040 ms0.5043 ms0.14
FloydWarshall_031200185.628 ms2.1240 ms1.9868 ms0.11
FloydWarshall_00240014,533.335 ms63.3518 ms52.9016 ms1.00
FloydWarshall_01 24006,507.878 ms28.2317 ms26.4079 ms0.45
FloydWarshall_02 24002,076.403 ms20.8320 ms19.4863 ms0.14
FloydWarshall_0324002,568.676 ms31.7359 ms29.6858 ms0.18
FloydWarshall_004800119,768.219 ms181.5514 ms160.9406 ms1.00
FloydWarshall_01 480055,590.374 ms414.6051 ms387.8218 ms0.46
FloydWarshall_02480015,614.506 ms115.6996 ms102.5647 ms0.13
FloydWarshall_03480018,257.991 ms84.5978 ms79.1329 ms0.15

Заключение

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

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

Надеюсь, вам было интересно посмотреть, как можно «поработать» над алгоритмом, в котором всего-то ТРИ цикла 🙂

Ссылки

  1. Floyd, R. W. Algorithm 97: Shortest Path / R. W. Floyd // Communications of the ACM. – 1962. – Vol. 5, №. 6. – P. 345-. https://dl.acm.org/doi/10.1145/367766.368168
  2. Ingerman, P. Z. Algorithm 141: Path matrix / P. Z. Ingerman // Communications of the ACM. – 1962. – Vol. 5, №. 11. – P. 556. — https://dl.acm.org/doi/10.1145/368996.369016

Автор материала – .NET Engineer ISsoft Олег Карасик

2 КОММЕНТАРИИ

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

ОСТАВЬТЕ ОТВЕТ

Please enter your comment!
Please enter your name here