Ускорить огромный цикл (Корреляционные моменты)

[Vamireh]

Есть алгоритм, считает корреляционные моменты. Если кратко, то представляет собой:

Код:

    QVector<int> vector;
    const int size = vector.size();

    double middle = 0.0;
    
    const auto end = vector.constEnd();
    for (auto it = vector.constBegin(); it != end; ++it)
        middle += *it;
    
    middle /= vector.size();
    
    for (int tau = tauMin; tau < tauMax; tau += step) {
        double y = 0;
        const auto iEnd = end - tau;

        for (auto it = vector.constBegin(); it != iEnd; ++it)
            y += ((*it) - middle) * ((*(it + tau)) - middle);

            y /= (size - tau - 1);
            result.append(QPointF(tau, y));
    }

Проблема в том, что размер vector может быть 4-5 млн., а tau от 0 до 100-200 тыс. В итоге через полчаса работы я вырубил программу. Можно ли как-то ускорить этот алгоритм до вменяемых времен расчета?

[Maestro]

Навскидку, попробуйте в начале циклов для result сделать reserve(int size);

[Igors]

Краткий и банальный ответ "никак". А чисто практически

1) *(it + tau) вылазит за границы массива

2) предраспределить result. Это не решит проблему но нужно для 3

3) очевидно распараллелить. С OpenMP дело сводится к добавлению одной строки

4) Если такие расчеты неизбежны - нужно делать их адаптивно и отображать в UI - напр через каждые 1000 точек. Еще лучше считать не подряд, а "решеткой". При наличии обновляемого графика юзверь занят наблюдением и время для него течет куда быстрее. Хороший способ "превратить отходы в доходы"

[Akon]

Действительно, *(it + tau) вылазит за границы массива.

1. Распараллелить, OpenMP как нельзя к месту.

2. Использовать SIMD инструкции современных процессоров.

Если допустимо, double заменить на float, - это даст ускорение в два раза. Если допустимо, перейти к арифметике с фиксированной точкой на int16, - это даст ускорение еще в два раза (относительно float).

С-либы для работы с SIMD - Intel IPP/MKL, у AMD - свои. Можно и на асме/интринсиках сделать, тем более, что алгоритм простой.

Чтобы вам лучше представить, что даст SIMD для вашей задачи, приведу пример: AVX инструкция DPPS оперирует 2-мя исходными 256-битными операндами (например, в каждом по 8 float'ов) и производит попарное умножение, т.е. за раз получается 8 вот таких результатов (правда без сложения, но это мелочь):

Код:

  // это один результат
        for (auto it = vector.constBegin(); it != iEnd; ++it)
            y += ((*it) - middle) * ((*(it + tau)) - middle);

Ввиду того, что АЛУ конвейеризовано, и имеется несколько исполнительных портов, то будет еще ускорение при последовательном исполнении однотипных инструкций (с независимыми операндами), например, упомянутая выше одна dpps занимает 12 тактов, а две - только 14. Ваш алгоритм хорош для этой возможности.

3. GPGPU (OpenCL, CUDA, ...). Ввиду того, что объем данных большой, возможно, будет лучше SIMD. Но здесь я не специалист, пока.

[Igors]

Если допустимо, double заменить на float, - это даст ускорение в два раза.

Моя практика этого не подтверждает, скорость остается примерно одинаковой. Возможно Вы имели ввиду "меньше объем данных для чтения".

Если допустимо, перейти к арифметике с фиксированной точкой на int16, - это даст ускорение еще в два раза (относительно float).

В первый раз слышу о такой - хранить флоты 16 бит (half) да, но об арифметике мне ничего не известно. Сообщите что/откуда. Спасибо

Чтобы вам лучше представить, что даст SIMD для вашей задачи, приведу пример: AVX инструкция DPPS оперирует 2-мя исходными 256-битными операндами (например, в каждом по 8 float'ов) и производит попарное умножение, т.е. за раз получается 8 вот таких результатов (правда без сложения, но это мелочь):

Имею эту радость в одном из проектов. Да, прирост скорости впечатляет. Но Вы забыли о др стороне медали - никто не обещал наличие AVХ на целевой машине. Там может быть по минимуму SSE3 а может и вообще ничего (хоть и редко). В результате приходится иметь 4 либы (для всех вариантов) и притыкивать их при старте. Не смертельно, но кол-во забот нарастает (особенно c учетом кросс-платформ)

[Old]

В первый раз слышу о такой - хранить флоты 16 бит (half) да, но об арифметике мне ничего не известно. Сообщите что/откуда. Спасибо

Вы занимаетесь 3d и никогда не слышали про арифметику с фиксированной запятой?  Впечатляет. Для реалтайма в далекие времена использовали его во все поля.
http://ru.wikipedia.org/wiki/%D0%A7%D0%B8%D1%81%D0%BB%D0%BE_%D1%81_%D1%84%D0%B8%D0%BA%D1%81%D0%B8%D1%80%D0%BE%D0%B2%D0%B0%D0%BD%D0%BD%D0%BE%D0%B9_%D0%B7%D0%B0%D0%BF%D1%8F%D1%82%D0%BE%D0%B9

Но врядли это здесь поможет, только если использовать 64битные числа для хранения.

[Akon]

Моя практика этого не подтверждает, скорость остается примерно одинаковой. Возможно Вы имели ввиду "меньше объем данных для чтения".

Ускорение за счет большего количества опреандов на регистр. Кол-во операндов обратно пропорционально размеру операнда. т.е. 4 double = 8 float. 

Имею эту радость в одном из проектов. Да, прирост скорости впечатляет. Но Вы забыли о др стороне медали - никто не обещал наличие AVХ на целевой машине. Там может быть по минимуму SSE3 а может и вообще ничего (хоть и редко). В результате приходится иметь 4 либы (для всех вариантов) и притыкивать их при старте. Не смертельно, но кол-во забот нарастает (особенно c учетом кросс-платформ)

Да, так и есть. Например, AVX начинается только с SandyBridge. Упомянутые либы от Intel имеют набор всех реализаций и выбирают оптимальный в рантайме.

[Igors]

Вы занимаетесь 3d и никогда не слышали про арифметику с фиксированной запятой?  Впечатляет. Для реалтайма в далекие времена использовали его во все поля.
http://ru.wikipedia.org/wiki/%D0%A7%D0%B8%D1%81%D0%BB%D0%BE_%D1%81_%D1%84%D0%B8%D0%BA%D1%81%D0%B8%D1%80%D0%BE%D0%B2%D0%B0%D0%BD%D0%BD%D0%BE%D0%B9_%D0%B7%D0%B0%D0%BF%D1%8F%D1%82%D0%BE%D0%B9

Но врядли это здесь поможет, только если использовать 64битные числа для хранения.

Выглядит как-то наивно - что я сохраню в 16 битах, по 1 цифре на целую/дробную?  А хоть бы и 64 - как тогда их складывать, "эмулировать" что ли 

Если уж разговор от том "кто чего слышал" - то вот есть формат 4 байта на точку (32-бита), который хранит RGB как 3 float значения (может быть напр тысячи, миллионы). Довольно популярен. Ну как он называется - можете и не знать, а вот как он хранит - догадаетесь? 

[Igors]

Ускорение за счет большего количества опреандов на регистр. Кол-во операндов обратно пропорционально размеру операнда. т.е. 4 double = 8 float. 

Уже больше 10 лет ассемблером не занимался (к сожалению). И продолжаю по старинке думать что регистры FPU 64-битные. Команда загрузки float из памяти в регистр сама  конвертирует  в double и наоборот. Или это уже не так? Спасибо

[Igors]

Да, так и есть. Например, AVX начинается только с SandyBridge. Упомянутые либы от Intel имеют набор всех реализаций и выбирают оптимальный в рантайме.

Хммм... сомневаюсь. Я использую embree от того же Intel и ничего она не выбирает - под что откомпилил - то и будет. На всякий случай спросил у авторов - да, компилите под каждый набор. А если "выбирает" то интересно как она это делает. Напр

Код:

с = a + b;   // хотя бы

Переменные должны быть обычными float или блоками по 4 float или даже по 8. И инструкции сложения тоже разные. И что, откомпилятся все варианты и в runtime код прыгнет на нужную ветку?  Хорошо бы конечно, но слабо верится, более реалистично есть какая-то директива указывающая компилятору под что делать код, ну или переменная окружения. Но "все в одном" так не выходит.

[Old]

Выглядит как-то наивно - что я сохраню в 16 битах, по 1 цифре на целую/дробную?   А хоть бы и 64 - как тогда их складывать, "эмулировать" что ли  

http://habrahabr.ru/post/131171/

как 3 float значения (может быть напр тысячи, миллионы).

Что имеется ввиду: каждая компонента может хранить число от 0 и до "тысячи, миллионы"?

[Akon]

Уже больше 10 лет ассемблером не занимался (к сожалению). И продолжаю по старинке думать что регистры FPU 64-битные. Команда загрузки float из памяти в регистр сама  конвертирует  в double и наоборот. Или это уже не так? Спасибо

Регистры FPU 80-битные и соответствуют типу long double. У SIMD расширений свое АЛУ (и свои регистры, за исключением MMX, которое использует регистры сопроцессора), к сопроцессору не имеющее никакого отношения.

Хммм... сомневаюсь. Я использую embree от того же Intel и ничего она не выбирает - под что откомпилил - то и будет. На всякий случай спросил у авторов - да, компилите под каждый набор. А если "выбирает" то интересно как она это делает.

Выбор в рантайме делается так:
- имеются реализации под все наборы инструкций в виде отдельных ДЛЛ;
- в рантайме на основании инструкции CPUID определяется набор поддерживаемых инструкций, ну и выбирается максимально функциональная (загружается соответствующая ДЛЛ).
Интерфейсные функции едины, просто внутри идет диспетчеризация в конкретную реализацию.

[Igors]

http://habrahabr.ru/post/131171/

Не думаю что сами Вы дочитали эту тягомотину до конца Да и к Вашей предыдущей ссылке это никакого отношения не имеет.

Что имеется ввиду: каждая компонента может хранить число от 0 и до "тысячи, миллионы"?

Да

По ходу дела заметим что прямолинейная попытка обрезать float (ну просто возьмем меньше бит на х-ку и мантиссу) ничего хорошего не приносит, half делает хитрее. Да, ну и все варианты скорость не повышают, а наоборот

[Igors]

Выбор в рантайме делается так:
- имеются реализации под все наборы инструкций в виде отдельных ДЛЛ;
- в рантайме на основании инструкции CPUID определяется набор поддерживаемых инструкций, ну и выбирается максимально функциональная (загружается соответствующая ДЛЛ).
Интерфейсные функции едины, просто внутри идет диспетчеризация в конкретную реализацию.

С этим я никак не спорю, т.к. сам делаю так же. Но я могу подгрузить достаточно большую/представительную единицу (напр ф-цию, класс), свободно писать арифметику и надеяться что "в runtime разберутся" возможности нет. Поэтому для задачи ТС все равно придется заводить неск компилируемых вариантов (что впрочем не смертельно)

[Old]

Не думаю что сами Вы дочитали эту тягомотину до конца

Конечно нет, я стараюсь не читаю хабр. Но это ваш уровень, там элементарные вещи на пальцах. Я бы порекомендовал faq фидо-эхи demo.design, но боюсь что для вас это будет сложновато. Судя по перлу ниже...

Да и к Вашей предыдущей ссылке это никакого отношения не имеет.

Да вы что правда? Вы и читать уже не умеете? Или это сложная тема для вас? Брысь учить арифметику и не позориться.