Какова вычислительная стоимость

Одной из основных проблем, с которой нам приходится сталкиваться при молекулярном моделировании, является расчет сил, зависящих от расстояния. Если мы можем ограничить функции силы и расстояния четными степенями расстояния разделения , то мы можем просто вычислить квадрат расстояния и не беспокоиться о . Однако, если существуют нечетные степени, нам нужно иметь дело с .r 2 = r ⋅ r r r = $r$$r^2 = {\bf r \cdot r}$$r$$r = \sqrt{r^2}$

Мой вопрос: насколько дорого обходятся вычисления реализованные в библиотеках общих языков (C / C ++, Fortran, Python) и т. Д.? Неужели действительно нужно много улучшений производительности при ручной настройке кода для конкретных архитектур?$\sqrt{x}$

В качестве дополнения к ответу moyner в , то на чипе sqrtобычно является rsqrt, то есть ответная квадратный корень , который вычисляет → 1 / √ . Так что, если в вашем коде вы собираетесь использовать только1/r(если вы занимаетесь молекулярной динамикой, то это так и есть), вы можете вычислитьнапрямую и сэкономить себе деление. Причина, по которойвычисляется вместо,состоит в том, что итерация Ньютона не имеет делений, только сложения и умножения.$a \rightarrow 1/\sqrt{a}$$1/r$r = rsqrt(r2)rsqrtsqrt

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

Некоторые более современные архитектуры, такие как архитектуры IBM POWER , не предоставляют rsqrtсами по себе, но дают оценку с точностью до нескольких бит, например FRSQRTE . Когда пользователь звонит rsqrt, это генерирует оценку, а затем одну или две (столько, сколько требуется) итераций алгоритма Ньютона или Голдшмидта с использованием регулярных умножений и сложений. Преимущество этого подхода состоит в том, что этапы итерации могут конвейеризоваться и чередоваться с другими инструкциями без блокирования FPU (очень хороший обзор этой концепции, хотя и на более старых архитектурах, см. В диссертации доктора философии Рольфа Штребеля ).

Что касается потенциалов взаимодействия, sqrtоперации можно полностью избежать, используя полиномиальный интерполятор потенциальной функции, но моя собственная работа (реализованная в mdcore) в этой области показывает, что, по крайней мере на архитектурах типа x86, sqrtинструкция достаточно быстрая.

Обновить

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

В контексте моделирования молекулярной динамики или любого моделирования на основе частиц с взаимодействием, ограниченным срезом, многое можно извлечь из лучших алгоритмов поиска соседей. Если вы используете списки ячеек или что-то подобное, чтобы найти соседей или создать список Verlet , вы будете вычислять большое количество ложных парных расстояний. В наивном случае только 16% проверенных пар частиц будут фактически находиться на расстоянии отсечки друг от друга. Хотя для таких пар взаимодействие не рассчитывается, доступ к данным частиц и вычисление паразитного парного расстояния сопряжены с большими затратами.

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

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

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

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

Есть также некоторые Streaming SIMD инструкции , которые могут быть использованы на регистрах XMM для быстрых векторных вычислений , найденных здесь . Эти регистры довольно малы, но если у вас есть известное количество координат (скажем, трехмерная декартова система координат), они могут быть немного быстрее.

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

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

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

Проект GROMACS MD возник из идеи использовать детали формата IEEE с плавающей запятой для создания итерационной схемы Ньютона-Рафсона для вычисления приемлемого приближения к обратному квадратному корню (см. Приложение B.3 в http: / /www.gromacs.org/Documentation/Manual ), но нет используемых процессоров HPC, где GROMACS все еще использует эту идею.