НОВОСТИ [Перевод] Ускорение в 14 000 раз или Победа компьютерной науки

NewsBot
Оффлайн

NewsBot

.
.
Регистрация
21.07.20
Сообщения
40.408
Реакции
1
Репутация
0
Как разработчику научного ПО мне приходится много программировать. И большинство людей из других научных областей склонны думать, что программирование — это «просто» набросать код и запустить его. У меня хорошие рабочие отношения со многими коллегами, в том числе из других стран… Физика, климатология, биология и т. д. Но когда дело доходит до разработки ПО, то складывается отчётливое впечатление, что они думают: «Эй, что тут может быть сложного?! Мы просто записываем несколько инструкций о том, что должен сделать компьютер, нажимаем кнопку „Выполнить” и готово — получаем ответ!»

Проблема в том, что невероятно легко написать инструкции, которые означают не то, что вы думаете. Например, программа может совершенно не поддаваться интерпретации компьютером. Кроме того, , не выполнив её. И есть много, очень много способов сильно замедлить выполнение программы. В смысле… реально замедлить. Так замедлить, что выполнение займёт всю вашу жизнь или больше. Это чаще всего происходит с программами, которые написаны людьми без компьютерного образования, то есть учёными из других областей. Моя работа — исправлять такие программы.

Люди не понимают, что информатика учит вас теории вычислений, сложности алгоритмов, вычислимости (то есть можем ли мы действительно что-то вычислить? Слишком часто мы считаем само собой разумеющимся, что можем!) Информатика даёт знания, логику и методы анализа, помогающие написать код, который выполнится за минимальное количество времени или с минимальным использованием ресурсов.

Позвольте показать пример огромной оптимизации одного простого скрипта, написанного моим коллегой.

В климатологии мы много масштабируем. Мы берём показания температуры и осадков из крупномасштабной глобальной климатической модели и сопоставляем их с мелкомасштабной локальной сеткой. Допустим, глобальная сетка равна 50×25, а локальная — 1000×500. Для каждой ячейки сетки в локальной сетке мы хотим знать, какой ячейке сетки в глобальной сетке она соответствует.

Простой способ представить проблему — это минимизировать расстояние между L[n] и G[n]. Получается такой поиск:


для каждой локальной ячейки L:
для каждой глобальной ячейки G[j]:
вычислить расстояние между L и G[j]
найти минимальное расстояние в наборе L * G
вернуть индекс нимимума

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


для каждой локальной ячейки L: # Сделать L раз
для каждой глобальной ячейки G[j]: # Сделать L x G раз
вычислить расстояние между L и G[j] # Сделать L x G раз
найти минимальное расстояние в наборе d[i*j] # Прочитать G ячеек L раз (L x G)
найти индекс, для которого ячейка соответствует минимуму # Прочитать G ячеек L раз (L x G)

Код выглядит примерно так:


obs.lon = (obs.lons-dxy)) &
(gcm.lons.lats[,1] = (obs.lats-dxy)) &
(gcm.lons.lats[,2] code>

Похоже на простой алгоритм. «Просто» вычислить расстояния, а затем найти минимум. Но в такой реализации по мере роста числа локальных ячеек наша стоимость вычислений растёт на её произведение с числом ячеек глобальной сетки. Канадские данные ANUSPLIN содержат 1068×510 ячеек (в общей сложности 544 680). Предположим, что наш GCM содержит 50×25 ячеек (в общей сложности 1250). Таким образом, стоимость внутреннего цикла в «каких-то вычислительных единицах» равна:



b0624323188b89de2b8495e7fc72675b.svg


где члены
e1a229081e8db6ee98dfb79797b987dd.svg
— это константы, соответствующие стоимости вычисления расстояния между двумя точками, нахождения минимальной точки и нахождения индекса массива. На самом деле, мы не заботимся (особо) о постоянных членах, потому что они не зависят от размера входных данных. Так что можно просто сложить их и оценить стоимость вычисления;



e22378e730943a0a5adbee12a87804a1.svg


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

680 миллионов.

Кажется, много… хотя, кто его знает? Компьютеры ведь быстрые, верно? Если мы запустим наивную реализацию, в реальности она отработает чуть быстрее 1668 секунд, что составляет чуть меньше получаса.


> source('BCCA/naive.implementation.R')
500 1000 1500 2000 2500 3000 ... 543000 543500 544000 544500 [1] "Elapsed Time"
user system elapsed
1668.868 8.926 1681.728

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

Другая важная структура, которой мы не воспользовались — тот факт, что широты идут в измерении
817b92407f764f57af9226e50cc788fd.svg
, а долготы — в измерении
9b34c4da5c757d4982bbd1b6f2e8998a.svg
. Таким образом, вместо выполнения операции
58f559255775c1f08315dd0c7b3b1cc0.svg
раз мы можем выполнить её
a7302669d3b3a0caa6369acdc4a5fa20.svg
раз. Это огромная оптимизация.

Как это выглядит в псевдокоде?


Для каждого local[x]:
bisect_search(local[x], Global[x])

Для каждого local[y]:
bisect_search(local[y], Global[y])

возвращается 2d сетка результатов поиска для каждого измерения

В коде:


## Perform a binary search on the *sorted* vector v
## Return the array index of the element closest to x
find.nearest 2)
if (x == v[mid]) {
return(mid)
} else if (x < v[mid]) {
return(find.nearest(x, v[1:mid]))
}
else {
return((mid - 1) + find.nearest(x, v[mid:length(v)]))
}
}

regrid.one.dim code>

Стоимость каждого поиска методом деления пополам — это log от входного размера. На этот раз наш входной размер разделён на пространство X и Y, поэтому будем использовать
5c8754b4832bde75f5580251bd545f7f.svg
и
06e67ded5cb0864af65332f98456f60c.svg
для Global, Local, X и Y.



a588fbfaeb8a091918960f99416e43c0.svg


Затраты получаются в 553 076. Что ж, 553 тысячи звучит намного лучше, чем 680 миллионов. Как это отразится на времени выполнения?


> ptm str(rv)
List of 2
$ xi: num [1:1068, 1:510] 15 15 15 15 15 15 15 15 15 15 ...
$ yi: num [1:1068, 1:510] 13 13 13 13 13 13 13 13 13 13 ...
>

0,117 секунды. То, что раньше выполнялось почти полчаса, теперь занимает чуть больше одной десятой секунды.


> 1668.868 / .117
[1] 14263.83

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

Раньше скрипт работал настолько медленно, что сохранял результат на диск для проверки учёным вручную, прежде чем продолжить. Теперь всё вычисляется в мгновение ока. Такие расчёты мы проводим сотни раз, так что в итоге экономим дни и даже недели вычислительного времени. И получаем возможность лучше взаимодействовать с системой, так что рабочее время учёных проходит с большей пользой… они не сидят без дела, ожидая окончания вычислений. Всё готово сразу.

Должен подчеркнуть, что для этих эпических улучшений производительности не требуется покупка каких-то больших компьютерных систем, распараллеливание или увеличение сложности… на самом деле код более быстрого алгоритма даже более простой и универсальный! Эта полная победа достигнута просто благодаря вдумчивому чтению кода и наличию некоторых знаний о вычислительной сложности.
 
Сверху Снизу