Информация об изменениях

Сообщение 2D-Linq и оптимизация цифровых фильтров от 25.06.2018 11:16

Изменено 25.06.2018 16:28 Sinclair

Старый текст
Новый текст
Разница

2D-Linq и оптимизация цифровых фильтров

Этот пост, по-хорошему, заслуживает место во flame.comp/, т.к. побудительным мотивом изначально являлся вопрос alex_public, заданный в недрах флейма в том форуме:

Но могу для разнообразия подкинуть ещё один примерчик для понимания. Простейшая, но очень нужная во многих областях задачка: имеем двухмерный массив чисел и надо получить из него новый массив, в котором каждая точка является средним от своих соседей (4-ёх или 8-ми, не суть). Что об этом скажет linq? )))

(https://rsdn.org/forum/flame.comp/7137817.1)

Но, поскольку речь идёт об пределах возмжностей кода на C#, я всё же запощу это сюда — в интересах тех, кого тоже волнует вопрос о производительности "математического" кода на дотнете, а перечитывать топики глубиной за сотню постов терпения нет.

Итак, что об этом скажет Linq? Вопрос, понятное дело, с подковыркой — типа, да Linq же — это убогая надстройка над IEnumerable, которую потом и соплями прикрутили к RDBMS. А тут у нас — принципиально другой материал, неперечислимый и, в отличие от, скажем, графов, к перечислениям несводимый.

Лично моими устами об этом Linq скажет следующее:

int[,] FourNeighborAverage(int[,] data) => from d in data select (d[-1, 0] + d[1, 0] + d[0, 1] + d[0, -1]) / 4;

Как это работает?
Ну, для начала у нас есть public static класс Array2d, где определён вот такой вот extension метод:

public static T[,] Select<T>(this T[,] source, Kernel<T> kernel)

То есть выражение d => (d[-1, 0] + d[1, 0] + d[0, 1] + d[0, -1]) / 4 трактуется как Kernel<int>:

public delegate T Kernel<T>(IRelativeAccess2d<T> point);

В свою очередь, интерфейс IRelativeAccess2d<T> определён весьма незатейливо:

public interface IRelativeAccess2d<T>
{
   T this[int dx, int dy] { get; } 
}

Теперь всё становится понятным: в наш Kernel передаётся интерфейс для доступа к "соседним" относительно текущего элементам массива, что позволяет нам описывать выражения типа "среднее от своих четырёх соседей".

Реализовать наш метод Select не очень трудно. Наивная реализация выглядит так:

public static T[,] Select<T>(this T[,] source, Kernel<T> kernel)
{
  var h = source.GetLength(0);
  var w = source.GetLength(1);
  var result = new T[h, w];
  for (int i = 0; i < h; i++)
    for (int j = 0; j < w; j++)
      result[i, j] = kernel(new RelativeArrayAccess2d<T>(source, i, j));
  return result;
}

структура RelativeArrayAccess2d<T> слишком тривиальна, чтобы её приводить — она хранит ссылку на массив source и координаты текущего элемента; при обращении по this[dx, dy] мы прибавляем относительные координаты к текущим и возвращаем элемент.
Этого кода уже достаточно для того, чтобы реализовать копирование массива в стиле linq:

from d in data select d[0, 0];

К сожалению, первый же запуск с ядром усреднения вылетит с ArrayIndexOutOfBoundsException. Логично — мы пробуем обращаться за пределы исходного массива. Суровые программисты на плюсах на такие мелочи внимания не обращают — ну, или требуют от программиста ручного указания размеров "каёмки", в которой фильтр не работает.
Правильное решение этой задачи — довольно сложно: нам придётся учиться подпиливать фильтр на ходу таким образом, чтобы он корректно работал на границах массива (у элемента рядом со "стороной" прямоугольника не 4, а 3 соседа, а у углового элемента — только 2. Строго говоря, именно их и нужно усреднять).
Но пока что мы обойдёмся решением в стиле С++ — ограничением области вывода фильтра.
Сделать это не очень сложно — мы просто один раз вызовем kernel со специальной реализацией IRelativeAccess, которая запоминает, куда идёт обращение:

private sealed class KernelMeasure<T> : IRelativeAccess2d<T>
{
  ...
  T IRelativeAccess2d<T>.this[int x, int y]
  {
    get
    {
      xmin = Math.Min(xmin, x);
      xmax = Math.Max(xmax, x);
      ymin = Math.Min(ymin, y);
      ymax = Math.Max(ymax, y);
      return default(T);
    }
  }
}

Вообще говоря, это будет работать не для всех типов Kernel — например, мы могли использовать Math.Random() для выбора смещения. Но для практически интересных случаев, когда смещения всегда постоянны, этот подход сработает.
Итого, у нас получается более-менее работоспособный код, который позволяет удобно записывать формулы фильтрации:

public static T[,] Select<T>(this T[,] source, Kernel<T> kernel)
{
  var h = source.GetLength(0);
  var w = source.GetLength(1);
  var km = new KernelMeasure<T>();
  kernel(km); // провели измерения
  var t = 0 - km.xmin;
  var b = h - km.xmax;
  var l = 0 - km.ymin;
  var r = w - km.ymax;
  var result = new T[h, w];

  for (int i = t; i < b; i++)
    for (int j = l; j < r; j++)
      result[i, j] = kernel(new RelativeArrayAccess2d<T>(source, i, j));
  return result;
}

На этом можно было бы и остановиться, если бы не ~~ехидный~~ комментарий в том же флейме на 30 уровней выше:

Ну как же, если мы говорим про работу с коллекциями (а не про СУБД), то обычный императивный код (типа for() ... ) никуда не делся. И при этом он по прежнему является как наиболее обобщённым, так и наиболее производительным.

CHALLENGE ACCEPTED!
По поводу "обобщённости" спорить смысла нет — очевидно, что приведённый from d in data select... можно без изменений заставить работать не только с 2d-массивом, но и с произвольной структурой, которая представляет из себя двумерную сетку, индексированную целочисленными координатами. Например, мы можем работать с "изображением" из double с размерами миллион на миллион, то есть имея 10^12 точек или около 7 терабайт. Перебор их циклом for упрётся даже не во время, а в размер памяти и диска локальной машины. Так что мы пилим всё на тайлы, раскидываем по нескольким тысячам серверов, отправляем код для исполнения, дожидаемся, ~~мёрджим~~. Ну, то есть получаем опять хэндл чудовищного "холста" для дальнейших операций — например, edge-detect, а потом ещё свёртка, а потом понижение размерности, и уже в конце получаем что-нибудь приемлемое для памяти локальной машины. И всё — в компактном и удобном для чтения (и написания) виде.

Вопрос у нас к тому пункту, где про "наиболее производительным".
Давайте напишем тот идеальный императивный код, который бы написал на C# программист, не доверяющий всем этми новомодным линкам, с их "динамикой и рефлексией":

static int[,] PerfectFourNeighborAverage(this int[,] source)
{
  var h = source.GetLength(0);
  var w = source.GetLength(1);
  var result = new int[h, w];

  for (int i = 1; i < h-1; i++)
    for (int j = 1; j < w-1; j++)
      result[i, j] = (source[i, j-1] + source[i, j+1] + source[i-1, j] + source[i+1,j])/4;
  return result;
}

Возьмём случайный массив размером 1000*1000, и попробуем сравнить производительность нашего "обобщённого" кода со "специфичным":

BenchmarkDotNet=v0.10.14, OS=Windows 10.0.15063.1088 (1703/CreatorsUpdate/Redstone2)
Intel Core i7-6600U CPU 2.60GHz (Skylake), 1 CPU, 4 logical and 2 physical cores
Frequency=2742190 Hz, Resolution=364.6720 ns, Timer=TSC
  [Host]     : .NET Framework 4.6.1 (CLR 4.0.30319.42000), 32bit LegacyJIT-v4.7.2650.0
  DefaultJob : .NET Framework 4.6.1 (CLR 4.0.30319.42000), 32bit LegacyJIT-v4.7.2650.0

     Method |      Mean |     Error |    StdDev | Scaled | ScaledSD |
----------- |----------:|----------:|----------:|-------:|---------:|
    Perfect | 15.008 ms | 0.1706 ms | 0.1512 ms |   1.00 |     0.00 |
       Linq | 33.423 ms | 0.6599 ms | 0.8345 ms |   2.23 |     0.06 |

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

Поэтому в следующей серии я попробую показать, что можно с этим сделать, ничего не меняя в коде бенчмарка. То есть те самые преимущества, которые бесплатны для прикладного программиста — от него потребуется разве что перекомпилировать свой проект. Никаких хинтов в коде — ни императивного развёртывания циклов, ни ручной векторизации, ни даже умелой расстановки "декларативных" прагм.
Только хардкор, только голый from d in data select (d[-1, 0] + d[1, 0] + d[0, 1] + d[0, -1]) / 4.