Прошу помощи по обсчёту математики

Mishka> Копирование займёт время (и опять двойной вложенный цикл).
memcpy работает достаточно быстро.

Mishka> Пробег по этим двум векторам — ещё цикл.
Но один, вместо двойного - по матрице.

Mishka> Код, написанный ручками для каждого элемента — ни одного ветвления, столько же умножений и сложений с возможным уменьшеием слагаемых
Да, именно это. Только это препроцессор генерирует в *.с файл в виде ф-ии, потом на твой выбор - DLL или просто статическая библиотека.

Mishka> В исходный, как было одно время очень популярно на том FORTAN-е.
Нее - я ничего не говорил про то что in и out - это один аргумент ( вида ввод-вывод ). Я считал что их два различных. Так вот - out можно использовать в процессе вычисления. Для операций вида +=

yacc>> Ты меня не понял - я имею ввиду загрузку DLL с отображением ф-ии calc в память. А в ф-ии матрица самим набором вычислений представлена. Ее код генерирует отдельная программа-препроцессор.
Mishka> Это динамический вызов, который дороже статического.
Нет никакой разницы. После того как dll загружена - это просто вызов ф-ии. Неважно шла ли она сначала с запуском приложения или подгрузилась потом. Время тратиться только на загрузку и GetProcAddress - далее разницы нет.

Mishka> Вызовы. Циклы.
Цикл уже по векторам для получения результата. А вызов - да, не избежать.

AXT> А, в общем, я ступил, надо было сразу выложить. Умножаем столбец на матрицу.
Для первой.

29	55	74	84
74	74	0	-74
84	-29	-74	55
55	-84	74	-29

Реализация в лоб (не считая расходов на циклы)
R0 = 29*V0+55*V1+74*V2+84*V3;
R1 = 74*V0+74*V1+ 0*V2-74*V3;
R2 = 84*V0-29*V1-74*V2+55*V3;
R3 = 55*V0-84*V1+74*V2-29*V3;
16 умножений, 12 сложений.

Введём промежуточные переменные и немного соптимизируем
I0 = 29*V0; I1 = 55*V0;
J0 = 29*V0; J1 = 55*V1;
K0 = 74*V2;
L0 = 29*V3; L1 = 55*V3;

R0 = I0 + I1 + K0 + L0 + L1;
R1 = 74*(V0+V1-V3);
R2 = I0 + I1 - J0 - K0 + L1;
R3 = I1 - J0 - J1 + K0 - L0;
8 умножений, 14 сложений.

64	64	64	64
83	36	-36	-83
64	-64	-64	64
36	-83	83	-36

Реализация в лоб (не считая расходов на циклы)
R0 = 64*V0+64*V1+64*V2+64*V3;
R1 = 83*V0+36*V1-36*V2-83*V3;
R2 = 64*V0-64*V1-64*V2+64*V3;
R3 = 36*V0-83*V1+83*V2-36*V3;
16 умножений, 12 сложений.

Очевидные упрощения:
R0 = 64*(V0+V1+V2+V3);
R1 = 83*(V0-V3)+36*(V1-V2);
R2 = 64*(V0-V1-V2+V3);
R3 = 36*(V0-V3)-83*(V1-V2);

Можно ввести две промежуточные:
I0 = V0 - V3; I1 = V1 - V2;

R0 = 64*(V0+V1+V2+V3);
R1 = 83*I0+36*I1;
R2 = 64*(V0-V1-V2+V3);
R3 = 36*I0-83*I1;

6 умножений, 10 сложений

В принципе, если компилятор поддерживает преобразование дерева вычислений, то он это найдёт. В режиме C++, но в C стиле это будет выглядеть как:

1. inline void mul1( const int* v, int*& r )
2. {
3.     // no code would be generated, we just introduce aliases
4.     const int& V0 = v[ 0 ];
5.     const int& V1 = v[ 1 ];
6.     const int& V2 = v[ 2 ];
7.     const int& V3 = v[ 3 ];
8.  
9.     const int I0 = 29 * V0;
10.     const int I1 = 55 * V0;
11.     const int J0 = 29 * V0;
12.     const int J1 = 55 * V1;
13.     const int K0 = 74 * V2;
14.     const int L0 = 29 * V3;
15.     const int L1 = 55 * V3;
16.  
17.     r[ 0 ] = I0 + I1 + K0 + L0 + L1;
18.     r[ 1 ] = 74 * ( V0 + V1 - V3 );
19.     r[ 2 ] = I0 + I1 - J0 - K0 + L1;
20.     r[ 3 ] = I1 - J0 - J1 + K0 - L0;
21. }
22.  
23. inline void mul2( const int* v, int*& r )
24. {
25.     // no code would be generated, we just introduce aliases
26.     const int& V0 = v[ 0 ];
27.     const int& V1 = v[ 1 ];
28.     const int& V2 = v[ 2 ];
29.     const int& V3 = v[ 3 ];
30.  
31.     const int I0 = V0 - V3;
32.     const int I1 = V1 - V2;
33.  
34.     r[ 0 ] = 64 * ( V0 + V1 + V2 + V3 );
35.     r[ 1 ] = 83 * I0 + 36 * I1;
36.     r[ 2 ] = 64 * ( V0 - V1 - V2 + V3 );
37.     r[ 3 ] = 36 * I0 - 83 * I1;
38. }

В общем, это примерно то же самое, что и в указанном Татариным алгоритме Штрассена (уменьшение количества умножений), только там про матрицу ничего не известно, а здесь всё известно. Ну и малые размеры вектора и матрицы позволяют рекурсию заменить развёрткой. А inline вместо вызова ф-ции положит эти инструкции прямо в нужное место. Кода получится больше, но работать будет быстрее. Можно попросить компилятор не оформлять стек, тогда inline будет ещё быстрее.

PS Для мат преоброзований матрицы надо почитать стандарт и немного теории, т.к. это всё похоже на фильтры. А там нельзя просто так переставлять строки и столбцы, т.к. последствия могут быть весьма интересными.

yacc> memcpy работает достаточно быстро.
Это цикл с тем самым переходом. На 386 архитектуре он может быть выражен командой, которая занимает вовсе не один цикл. А копирование вектора N раз — это ещё один цикл и memcpy вызывается внутри этого цикла. Т.е. ты убрал цикл из умножение, но добавил в копирование.

yacc> Но один, вместо двойного - по матрице.

Он переместился наверх и вовсе не исчез.

yacc> Да, именно это. Только это препроцессор генерирует в *.с файл в виде ф-ии, потом на твой выбор - DLL или просто статическая библиотека.

Статическая лучше. Линкер подправит все адреса, они будут в водопроводе при чтении команды call. А в случае динамического вызова в водопроводе адрес, по которому лежит адрес, на который надо перейти. И предсказать очень трудно, т.к. команда перед call может поменять это значение.

yacc> Нее - я ничего не говорил про то что in и out - это один аргумент ( вида ввод-вывод ). Я считал что их два различных. Так вот - out можно использовать в процессе вычисления. Для операций вида +=

Зато я говорил. А отдельную память — всегда.

yacc> Нет никакой разницы. После того как dll загружена - это просто вызов ф-ии. Неважно шла ли она сначала с запуском приложения или подгрузилась потом. Время тратиться только на загрузку и GetProcAddress - далее разницы нет.

Блин, есть и большая. И поэтому виртуальные методы в том же C++ всегда медленней, чем обычные. GetProcAddress есть двух видов, вообще-то. Просто в одном случае таблица адресов инициализируется по номерам, а в другом случае по именам (текстовым константам). После этого все вызовыв выглядят так:

1.  
2. call  offset(table)

Где offset вычисляется компилятором исходя из номера процедуры, а table — адрес таблицы с количеством входов, описанным в заголовке DLL. Т.е. настоящий адрес процедуры лежит в табличке. А в случае статической библиотеки это дело имеет вид:

1.  
2. call  mul1

Где на момент стоят нули, в специально табличке в объектном файле есть вход с внешним именем (то, которое включает mangled name в случае C++) и смещение относительно кода, где надо подправить. Если вызовов 5, то будет 5 таких входов в таблички. В момент сбора программы, линкер поместит библиотеку по определённому адресу и всякий новый объектный файл будет обработан так, что все входы в табличку будут обработаны, все 0 будут замещены на адрес, где библиотека сидит (точнее процедура из библиотеки). Это всё в предположении, что программа будет загружена в память на исполнение с определённого базового адреса. Если базовый адрес изменится (в Виндах есть возможность программку сместить, а потом отобразить странички из других DLL или ядра в освободившееся пространство), то загрузчик пробежит по всем местам, где линкер бегал и поправит адреса — вся инфа есть в закогловке exe файла. Кстати, требования к COM были таковы, что код должен был быть relocatable, т.е. его можно было смещать в памяти без такого пробега. Поэтому можно было использовать команды только с относительной адресацией, а не с абсолютной, что ставило крест на длинных переходах. Для длинных переходов приходилось подсчитывать смещение на регистре или в памяти и использовать переход по регистру/памяти.

yacc> Цикл уже по векторам для получения результата. А вызов - да, не избежать.
Он есть до — копирование вектора N раз. Причём вложенный.

Mishka> Т.е. ты убрал цикл из умножение, но добавил в копирование.
Только он более быстрый - memcpy - оптимизированная.

Mishka> Статическая лучше. Линкер подправит все адреса, они будут в водопроводе при чтении команды call.
Так используй ее.

yacc>> Цикл уже по векторам для получения результата. А вызов - да, не избежать.
Mishka> Он есть до — копирование вектора N раз. Причём вложенный.
Ты не понял. Я здесь говорю уже о функции в которой уже нет никаких ветвлений. Вызов - для вызова ее - call mult - его не избежать. И цикл уже по векторам.

yacc> Только он более быстрый - memcpy - оптимизированная.
Она очень сильно оптимизированная в основном только на Интеле, бо есть аппаратная поддержка. У Microsoft-а, ЕМНИП, даже были с разной мнемоникой для байтов, слов и двойных слов.
yacc> Так используй ее.
DLL не позволяет по своей сути.

yacc> Ты не понял. Я здесь говорю уже о функции в которой уже нет никаких ветвлений. Вызов - для вызова ее - call mult - его не избежать. И цикл уже по векторам.
Блин, цикл-то сам по себе никуда не делся и будет выполняться. Внутри ф-ции или строго перед её вызовом — уже и не важно, т.к. каждый вызов влечёт такую подготовку.

А вызов избежать можно при помощи inline.

Вот, что надо сделать, это посмотреть, как всякие SIMD инструкции могут повлиять. Или DSP. Особенно, если ядер много. Из серии, все умножения вектора на матрицу (16) можно сделать сделать 4 устройствами за один такт, а результат положить так, чтобы за такт сделать сложения сделать за 2 такта (все V0+V1 на одном, V2+V3 на другом за такт), а потом сложить суммы ещё раз на одном устройстве за такт. Но надо считать сколько заёмёт подготовка.

Mishka> Я работал с DSP очень мало, но цикл и там не самая дешёвая вещь. Даже аппаратно поддержанный. Выпрямленный цикл, если память позволяет всё вместить, быстрее.

It depends. На том же BlackFin аппаратный цикл стоит 1 такт на запуск (или два всё же? подзабыл), остальные итерации бесплатны, плюс 2 вложенных цикла аппаратно поддерживаются. Но надо следить за конвейером, а то будет дико тормозить на самом деле.

Mishka> В принципе, если компилятор поддерживает преобразование дерева вычислений, то он это найдёт. В режиме C++, но в C стиле это будет выглядеть как:

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

Mishka> А inline вместо вызова ф-ции положит эти инструкции прямо в нужное место. Кода получится больше, но работать будет быстрее. Можно попросить компилятор не оформлять стек, тогда inline будет ещё быстрее.

Стека нет вообще, ни возвратов, ни данных Это ж числодробилка, нафига он ей?

Mishka> PS Для мат преоброзований матрицы надо почитать стандарт и немного теории, т.к. это всё похоже на фильтры. А там нельзя просто так переставлять строки и столбцы, т.к. последствия могут быть весьма интересными.

Тут, внезапно, счастье: математика по стандарту целочисленная. А правила округления при делении — стандартные, выкидываем дробную часть. Халява, сэр

yacc>> Только он более быстрый - memcpy - оптимизированная.
Mishka> Она очень сильно оптимизированная в основном только на Интеле, бо есть аппаратная поддержка.

... причём первым её сделал, как ни странно, AMD в своём Athlon Там при счётчике более ЕМНИМС 64 включался аппаратный копировщик.

Mishka> У Microsoft-а, ЕМНИП, даже были с разной мнемоникой для байтов, слов и двойных слов.

memcpy или команда процессора? Команды то всю дорогу были с указанием размера rep movsb, rep movsw, rep movsd

Mishka> Посмотрел Инетеловские доки. Умножение стало почти так же быстро работать на последних сериях i7, как и сложение. А вот деление ...

(вспомнилась что-то история про Pentium divide bug, но лень рассказывать)

Самое быстрое деление — по табличке Ну или начальное по таблице, а потом Ньютоном добить, как в 3DNow!/SSE/AVX. Сам грешен, делал такое в железе
Кстати, в случае с плавающей запятой есть интересный момент: так как числа нормализованы, то деление на 2^x * Q, где Q — малое целое число, сюрприз, не требуют коррекции Ньютоном вообще — сразу выходит точное значение.

Mishka> Она очень сильно оптимизированная в основном только на Интеле, бо есть аппаратная поддержка.
Кто мешает оптимальный алгоритм сделать для другого процессора?

yacc>> Так используй ее.
Mishka> DLL не позволяет по своей сути.
Нет - статическую линковку. Сгененрируй файл, скомпилируй и подлинкуй непосредственно к приложению, без оформления в виде DLL.

Mishka> А вызов избежать можно при помощи inline.
Когда данные читаются все равно цикл остается. Инлайн тебе не даст гибкости если несколько константных матриц используются - т.е. есть ф-ии calc1 ... calcN.

Mishka> Вот, что надо сделать, это посмотреть, как всякие SIMD инструкции могут повлиять. Или DSP. Особенно, если ядер много.
Или регистры общего назначения - там тоже можно оптимизировать.

yacc> Когда данные читаются все равно цикл остается. Инлайн тебе не даст гибкости если несколько константных матриц используются - т.е. есть ф-ии calc1 ... calcN.
Если внутри ф-ции, то это неверно.