Дмитpий Hecтepук

Блог о программировании — C#, F#, C++, архитектура, и многое другое

Снова про оптимизацию математических выражений

11 комментариев

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

Инлайнинг возведения в степень

Если ваша программа реализует f(x)=ax^2+bx+c как

double f(double a, double b, double c, double x)
{
  return a * std::pow(x, 2.0) + b*x + c;
}

то вы оказываете пользователю медвежью услугу. Для начала, следует помнить, что возведение в целочисленную степень является достаточно бессмысленной операцией, т.к. простое умножение “инлайн” или умножение в цикле всегда будет быстрее. Причина этого проста – функции вроде std::pow() расчитаны на подсчет степени для любых показателей, в т.ч. нецелочисленных. Как только в экспоненту ставится число вроде 2.0, строятся ряды и производительность идет лесом.

Решение этой проблемы весьма простое. Во-первых, для всех простых показателей вроде x^2, x^3 и так далее лучше все инлайнить, т.е. писать x*x, x*x*x и т.д. Что касается более крупных степеней то тут тоже все просто – надо чтобы степень была целым числом, т.е. имела тип int. Тогда будет вызвана пегрузка std::pow(), которая в свою очередь вызовет функцию _Pow_int() (по крайней мере в <math.h> от Microsoft), и все будет хорошо.

В .NET все не так ажурно т.к. Math.Pow() берет два параметра типа double. Очень досадное упущение, но это можно устранить самостоятельно:

public static long IntPower(int x, short power)
{
  if (power == 0) return 1;
  if (power == 1) return x;
  int n = 15;
  while ((power <<= 1) >= 0) n--;
  long tmp = x;
  while (--n > 0)
    tmp = tmp * tmp *
         (((power <<= 1) < 0) ? x : 1);
  return tmp;
}

Пример выше – для типа short, т.к. его обычно вполне достаточно. Естественно, следует также понимать что есть некий порог, после которого дешевле использовать канонический pow(double,double) вместо итеративного умножения.

Сокращение количества умножений

Допустим что мы воспользовались советом выше, и наша функция теперь выглядит вот так:

double f(double a, double b, double c, double x)
{
  return a*x*x + b*x + c;
}

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

double f(double a, double b, double c, double c)
{
  return x*(a*x + b) + c;
}

На самом деле, факторизация выражений, т.е. вынесение за скобки общих членов — задача не из легких, т.к. требует чтобы система знала уйму всяких разных правил. Например, как записать f(x,y)=x^2-y^2? Вы действительно думаете, что без двух умножений не обойтись? Ничего подобного!

double f(double x, double y)
{
  return (x + y) * (x - y);
}

И если равенство a^2-b^2=(a+b)(a-b) мы еще можем запомнить, то какие-то более сложные равенства уже не получится, и придется жить без них. Хотя у пользователя всегда есть вариант использовать какой-нибудь полноценный математический пакет и собственноручно оптимизировать вычисление.

Предварительные вычисления

Понятное дело что выражение \frac{1}{3}x превратится во что-то вроде 0.333333333333 * x будучи оптимизированным компилятором. Но если взять выражение вроде \frac{2x}{3y}, то тут уже все неоднозначно. До какой степени компилятор сможет оптимизировать и “инлайнить” результаты вычислений? Я думаю что в репертуар типичного компилятора не входит глубокий анализ и оптимизация вычислений.

Выделение переменных

Одна из проблем с генерацией кода для вычислений (к пр. на основе Excel) заключается в том, что есть большой шанс выполнить одну и ту же дорогую операцию несколько раз. Например, если напрямую конвертировать выражение 2cos^3x-cos^2x, то в результате можно получить весьма предсказуемый и неэффективный код:

double f(double x)
{
  return 2.0*pow(cos(x),3.0) - pow(cos(x), 2.0);
}

А если еще “проинлайнить” степень, то будет вообще ужасно. Мораль в том, что при генерации подобных выражений, нужно уметь оценивать стоимость различных вычислений и не давать общим элементам с большой стоимостью повторяться. Например, код выше имеет смысл записать вот так:

double f(double x)
{
  double u = cos(x);
  return (2*u - 1)*u*u;
}

То, что приведено выше – это так называемая схема Горнера. Примечательно, что эта задача — анализ выражений на поиск общих кусков кода — будучи созданной ускорять код, сама по себе имеет огромную стоимость, т.к. подразумевает поиск любого поддерева с определенной стоимостью в остальных частях дерева и выделения их для последующей замены. В общем случае эта задача не решается. К счастью, при трансляции Excel→С++ ее можно решить хотя бы частично: результаты промежуточных вычислений часто хранятся в отдельных клетках, а детектировать их повторное использование легче т.к. это всего лишь сравнение cell reference’ов.

Заключение

Для создания вменяемых математических моделей недостаточно просто конвертировать код “один в один”. В конце концов, главная цель подобных конверсий из Excel в С++ — получить многократный прирост производительности. А сколько эта конверсия займет времени никого толком не волнует.

Реклама

Written by Dmitri

12 февраля 2012 в 22:26

Опубликовано в Mathematics

комментариев 11

Subscribe to comments with RSS.

  1. Гляньте вот в эту статью http://www.codeproject.com/Articles/87294/Symbolic-Calculation-in-F
    Думаю, подобным образом можно обучить код делать препроцессинг. Конечно, если сами формулы hard-coded.

    Vlad

    13 февраля 2012 at 2:18

    • Ну да, у меня в MathSharp как раз построен такой большой discriminated union, который покрывает все поддерживаемые выражения. К сожалению, модель с которой я работаю в Excel уже написана через API на C# — он использует обычный полиморфизм и с ним намного сложнее работать.

      Dmitri

      13 февраля 2012 at 10:15

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

        А вот если она (формула) записывается в исходном коде, то вполне такая интересная задача получается. На автоматическое упрощение выражения.
        Только я не очень представляю, как реализовывать, например, метод добавления и вычитания (или домножения и деления) одинаковых слагаемых.
        Или выделение полного квадрата.
        Или поиск производной, и если она =0, то само выражение — всегда константа.
        И всё такое.

        Vlad

        13 февраля 2012 at 11:18

        • На самом деле, следует помнить что и в ReSharper уже есть некоторые механизмы по упрощению вычислений. Например, у меня недавно MathSharp сгенерировал что-то вроде 1.0/2.0, и R# предложил мне контекстное действие которое заменило это на 0.5. Естественно, это лишнее телодвижение, но всяко лучше чем ничего.

          Dmitri

          14 февраля 2012 at 13:58

        • Вот когда они научатся понимать myVar-myVar — тогда от этого будет польза, да.
          Следующий шаг — научиться понимать Math.Pow(Math.Sin(myVar), 2.0) + Math.Pow(Math.Cos(myVar), 2.0)
          А константы заменять — согласитесь, не фича. :)

          Vlad

          14 февраля 2012 at 15:30

        • Как раз это сейчас можно хоть как-то понять… конечно не как sin^2x+cos^2x=1, а скорее как тот факт что Math.Pow(Z,2.0) == Z*Z ну и соответственно то что нужно ввести переменные, т.е.

          var s = Math.Sin(myVar);
          var c = Math.Cos(myVar);
          return s*s + c*c;
          

          Конечно это не так круто как просто написать return 1, но мы пока до такого не дотягиваем.

          Dmitri

          14 февраля 2012 at 15:56

        • Я про MathSharp и говорил.
          Если в общем случае (для любого кода) это сделать слишком сложно, то Вам проще: у Вас есть на входе одно выражение, и в рамках этого выражения значение переменной не меняется.

          И потому Вы можете запилить паттерны, с помощью которых можно упростить выражение.

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

          Может, Вам, кстати, будет интересен ещё один человек http://nogin.org/ и его проект http://metaprl.org/default.html
          Мне он, правда, так и не ответил, но разработка интересная. Она умеет делать математические пруфы.

          Vlad

          14 февраля 2012 at 16:32

  2. Дмитрий, а все эти вещи как раз прекрасно дают в ВУЗе (я про математику ;)) Это к спору о необходимости получения ВО

    Рашид Фатыхов

    13 февраля 2012 at 7:08

    • Нет, их дают в школе. До ВУЗа многое уже забывается.

      Dmitri

      13 февраля 2012 at 10:08

      • Ага, и ряды дают в школе? )))

        Рашид Фатыхов

        13 февраля 2012 at 10:31

        • Насчет рядов не знаю, может и не дают… но мне и в университете не давали т.к. я учился на comp sci а не на математике.

          Dmitri

          14 февраля 2012 at 13:55


Оставить комментарий

Заполните поля или щелкните по значку, чтобы оставить свой комментарий:

Логотип WordPress.com

Для комментария используется ваша учётная запись WordPress.com. Выход / Изменить )

Фотография Twitter

Для комментария используется ваша учётная запись Twitter. Выход / Изменить )

Фотография Facebook

Для комментария используется ваша учётная запись Facebook. Выход / Изменить )

Google+ photo

Для комментария используется ваша учётная запись Google+. Выход / Изменить )

Connecting to %s

%d такие блоггеры, как: