[image]

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

Особенно к Мишке и 101 просьба
 

AXT

инженер вольнодумец
★★
Преамбула:

Надо считать математику. У меня болезнь и соответственно тупняк, так что прошу напомнить куда копать. Ибо работать надо. Спецов лучше меня в зоне досягаемости всё равно нет, так что решать математику придётся мне.

То есть, решалка в принципе есть, но ползает со скоростью умираюхей черепахи. Для исполнения ТЗ требуется ускорение примерно в 30 раз.

Амбула:

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

И тут грабли. Стандартный фокус с формированием матрицы из векторов и упрощённым умножением не катит, т.к. матрица всего 32x32. Толку нет (хотя у меня есть совсеи уж отчаянные мысли, но там вообще жесть будет).

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

Кто что посоветует?
   
RU Massaraksh #12.12.2014 00:41  @Sandro#10.05.2014 21:31
+
-
edit
 

Massaraksh

аксакал
★☆
AXT> Если выкинуть всякую малозначащую труху, то мне надо умножать вектор на матрицу. Много раз.

Миллион раз в секунду?
   33.033.0
RU Серокой #12.12.2014 01:00  @Sandro#10.05.2014 21:31
+
-
edit
 

Серокой

координатор
★★★★
Тупо: а мультимедийные расширения процессора нельзя дёрнуть? SSE?
   

yacc

старожил
★★★
AXT> Кто что посоветует?
Попробуй использовать LAPACK - там вроде должны быть алгоритмы оптимизированные.
Библиотека есть для Си и Фортрана точно.
   39.0.2171.7139.0.2171.71
RU AXT #12.12.2014 06:08  @Серокой#12.12.2014 01:00
+
-
edit
 

AXT

инженер вольнодумец
★★
Серокой> Тупо: а мультимедийные расширения процессора нельзя дёрнуть? SSE?

Да я как раз и думаю в эту сторону, но есть non-intel платформы, плюс (самое страшное) данные лежат не подряд.
   13.0.782.22013.0.782.220
+
-
edit
 

AXT

инженер вольнодумец
★★
Massaraksh> Миллион раз в секунду?

Около 20 млн в идеале. Матрица 4*4 — немного, но надо же синхроизироваться, и вот тут залежи граблей ...

(правка) 4*4 она после беспощандной редукции... похоже дальше никак.
   13.0.782.22013.0.782.220
RU Massaraksh #12.12.2014 23:37  @Sandro#12.12.2014 06:11
+
-
edit
 

Massaraksh

аксакал
★☆
Massaraksh>> Миллион раз в секунду?
AXT> Около 20 млн в идеале. Матрица 4*4 — немного, но надо же синхроизироваться, и вот тут залежи граблей ...
AXT> (правка) 4*4 она после беспощандной редукции... похоже дальше никак.

Крайгинг-интерполяция?
А если на матрицу сделать маску, и не проводить лишних умножений?
   
EU Татарин #15.12.2014 13:57  @Sandro#10.05.2014 21:31
+
-
edit
 

Татарин

координатор
★★★★★
AXT> Но матрица имеет совсем немного значений элементов и почти регулярную структуру, так что есть подозрения, что там скрыт более простой алгоритм. Ранг точно меньше 32, и сильно.
AXT> Кто что посоветует?
Первое что сделать - выделить линейно-независимые элементы, поискать в матрице какие-то симметрии, пусть неочевидные (симметрия с перестановками столбцов/строк - тоже симметрия). Ну и т.п. - поиграться с математикой, короче.
Можно даже втупую: выписать (машиной, конечно :)) полные выражения для элементов результирующего тензора, подставив матрицу как константы и упростить их (поэлементно) машиной.

На математике можно выиграть порядки.
А на библиотеках/SIMD - в лучшем случае только разы.

Если математика не помогает, то я б смотрел на CUDA/OpenCL (в зависимости от железа). Там такие задачи с массивным перемножением всего на всё бегут замечательно.
   39.0.2171.9539.0.2171.95
Это сообщение редактировалось 15.12.2014 в 14:05
US Mishka #16.12.2014 22:18  @Татарин#15.12.2014 13:57
+
-
edit
 

Mishka

модератор
★★★
AXT,

Пардон, не видел. Татарин уже ответил по существу. Я только его ответ прокомментирую.

Татарин> Первое что сделать - выделить линейно-независимые элементы, поискать в матрице какие-то симметрии, пусть неочевидные (симметрия с перестановками столбцов/строк - тоже симметрия). Ну и т.п. - поиграться с математикой, короче.

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

Татарин> Можно даже втупую: выписать (машиной, конечно :)) полные выражения для элементов результирующего тензора, подставив матрицу как константы и упростить их (поэлементно) машиной.

Даже матрица 32 на 32 и вектор на 32 — уже можно сделать 32 линейные записи по 32 элемента в каждом — распрямить циклы. Если в матрице много 0, 1, 2, то эти выражения вместо умножения легко оптимизируются ручками (для 2 в целочисленной арифметике). Если не делал преобразований матрицы, а есть одинаковые строчки (для вектора столбца), то имеет смысл сравнить значения соответствующие в векторе и не перевычислять.

Татарин> На математике можно выиграть порядки.

Я бы посоветовал погонять умножения и разные запросы в Maple. Она может предоставить совершенно нетривиальный вариант упрощения умножения для конкретной матрицы (матрица, как я понял, константная).




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

AXT, а какие данные (целочисленные? с плавающей? с фиксированной? диапазон значений?), какая точность (если с фиксированной, то может легче смоделировать целочисленной?)?
   34.034.0
EE Татарин #17.12.2014 02:45  @Mishka#16.12.2014 22:18
+
-
edit
 

Татарин

координатор
★★★★★
Mishka> AXT, а какие данные (целочисленные? с плавающей? с фиксированной? диапазон значений?), какая точность (если с фиксированной, то может легче смоделировать целочисленной?)?
К этому ещё добавлю - нельзя ли посмотреть на матрицу?
Если секретно-секретно, ну, сконструируй так похожее, чтоб суметь извлечь пользу от анализа коллективным разумом. Домножь там на коэффициент её, что ли. :D
   39.0.2171.9539.0.2171.95

yacc

старожил
★★★
Mishka> Даже матрица 32 на 32 и вектор на 32 — уже можно сделать 32 линейные записи по 32 элемента в каждом — распрямить циклы.
Там возможно проще N раз клонировать входной вектор в вектор с дубликатами. А константную матрицу также в виде вектора представить - NxN длиной. Потом просто один пробег попарного умножения плюс прибавления к результату с попутным вычислением целевой ячейки - чтобы она на N+1 элементе меняла значение.
Ну или просто попарное умножение, а потом N вызовов ф-ии сложения элементов вектора с указанием диапазона.
   39.0.2171.9539.0.2171.95
EE Татарин #17.12.2014 03:58  @yacc#17.12.2014 03:08
+
-
edit
 

Татарин

координатор
★★★★★
Mishka>> Даже матрица 32 на 32 и вектор на 32 — уже можно сделать 32 линейные записи по 32 элемента в каждом — распрямить циклы.
yacc> Там возможно проще N раз клонировать входной вектор
Для общего случая есть алгоритм Штрассена.
Но АХТ утверждает, что у него матрица какая-то неполноценная :) - с нулевыми элементами и т.п, и вообще константа. И тут на самом деле может быть, где гульнуть.
   39.0.2171.9539.0.2171.95

101

аксакал

Автор, кинь матрицу и вектора. Посмотрят наши, что нибудь посоветуют.
   1414

AXT

инженер вольнодумец
★★
101> Автор, кинь матрицу и вектора. Посмотрят наши, что нибудь посоветуют.

Сейчас уже не имеет смысла — ещё летом всё было заборото оптимизацией "в лоб", так что интересно только в смысле поиграться с математикой. Матрицы в стандарте ITU-T H.264, если нет доступа к, то могу выложить матрицы. Они там реально странные, вроде как закономерности есть даже на глаз, но единого принципа я так и не нашёл.

А сейчас попал на похожую задачу в рамках хобби :) Там тоже 4*4, но источник другой.

И всем, разумеется, спасибо за советы :hat-off:
   13.0.782.22013.0.782.220
RU yacc #17.12.2014 09:51  @Татарин#17.12.2014 03:58
+
-
edit
 

yacc

старожил
★★★
Татарин> И тут на самом деле может быть, где гульнуть.
Навскидку linear algebra - Recommendations for a usable, fast C++ matrix library? - Computational Science Stack Exchange
   39.0.2171.9539.0.2171.95
RU Massaraksh #17.12.2014 22:41  @Mishka#16.12.2014 22:18
+
-
edit
 

Massaraksh

аксакал
★☆
Mishka> Massaraksh, я думаю, что ручками оптимизоровать выражения при константной матрице будет куда более эффективно. При этом даже можно врубить для данноё ф-ции максимальную оптимизацию выражений при компилировании данной ф-ции по быстродействию в ущерб размеру. Она всё равно уложится в кэш у современных процессоров. Это будет лучше маски, которую надо прикладывать динамически.
Просто из контекста не совсем понятно, насколько она константная. У меня недавно была задача, как раз связанная с крайгинг-интерполяцией, где тоже приходилось много раз умножать вектор на константную матрицу, которая была константной аж 1/24 секунды.
   
Это сообщение редактировалось 18.12.2014 в 21:21

Mishka

модератор
★★★
yacc> Там возможно проще N раз клонировать входной вектор в вектор с дубликатами.

Зачем? Какой смысл? Если предпологать многопоточные процы, то векторы входной и на момент выполнения операцыии является тоже константным. По чтению нет никакой конкунренции. Или ты предлагаешь анализ вектора на дубликаты? Если последнее, то при разных строчках матрицы (для вертикального вектора, для горизонтального — стобцы матрицы), выход будет другим. 32 не настолко большая цифра, чтобы можно было бы выиграть на динамическом анализе, если вообще можно выиграть.

yacc> А константную матрицу также в виде вектора представить - NxN длиной. Потом просто один пробег попарного умножения плюс прибавления к результату с попутным вычислением целевой ячейки - чтобы она на N+1 элементе меняла значение.

Распрямение умножение — это 32 выражения вычисления, т.е. чисто представлено уже в виде таком. Если динамический анализ проведён, то я бы написал кучку разных алгоритмов для умножения именно на "такой" вектор, если он позволит сильно сократить количество умножений-сложыений (заметь, что в целофичленной арифметике сложение предопчительней умноженюю, а вот в плавающей ежу может быть наоборот — умножение не требует нормализации).

yacc> Ну или просто попарное умножение, а потом N вызовов ф-ии сложения элементов вектора с указанием диапазона.
Блин, зачем попарные умножения? Вызов стоит дорого. Ты представляешь, что такое динамический поиск общих подвыражений, которые ты пытаешься сделать?
   31.031.0
US Mishka #25.12.2014 13:45  @Massaraksh#17.12.2014 22:41
+
-
edit
 

Mishka

модератор
★★★
Massaraksh> Просто из контекста не совсем понятно, насколько она константная. У меня недавно была задача, как раз связанная с крайгинг-интерполяцией, где тоже приходилось много раз умножать вектор на константную матрицу, которая была константной аж 1/24 секунды.

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

yacc

старожил
★★★
yacc>> Там возможно проще N раз клонировать входной вектор в вектор с дубликатами.
Mishka> Зачем? Какой смысл?
Затем, что матрица последовательно лежит в памяти - дать ей такой же вектор для попарного умножения типа a[i] * b[i}

Mishka> Распрямение умножение — это 32 выражения вычисления
Когда я возился картами рельефа основная просадка была ... на ветвлении. Умножение выполняется быстро.

Mishka> Блин, зачем попарные умножения? Вызов стоит дорого. Ты представляешь, что такое динамический поиск общих подвыражений, которые ты пытаешься сделать?
Тогда самый быстрый способ - динамическая подлинковка.
Т.е. файл, который содержит ф-ию умножения на константную матрицу - отдельный.
Например mult.c содержит ф-ию f( value_type_ptr _inPtr, value_type_ptr _outPtr ) ( чисто к примеру )
А сам файл формируется динамически - т.е. некоторый препроцессор его формирует по константной матрице. сама ф-ия просто содержит набор вида _out[10] = _in[5]*10 + _in[6]*11 + ...
никаких ветвлений.
самое быстрое вычисление.
   

Mishka

модератор
★★★
yacc> Затем, что матрица последовательно лежит в памяти - дать ей такой же вектор для попарного умножения типа a[i] * b[i}

Ы? Что это даст? Какая выгода?

И последовательно лежит — это как? В том же FORTRAN-e матрица лежит по столбцам, а в PL/I, C, C++ — по строкам. В Алголе 68 она вообще может быть разбросана по памяти поэлементно. Ну и мы не знаем, что представляет собой вектор — столбец или строку. Но и это всё не ничего не даёт, т.к. вектор нет необходимости дублировать — пофигк умножению и сложению. И записывать обратно в вектор не можем до окончания операции полностью, т.к.о значение последнего элементо вектора-результата зависит от каждого элемента входного вектора.

yacc> Когда я возился картами рельефа основная просадка была ... на ветвлении. Умножение выполняется быстро.
Какое умножение? На каком процессоре? Мы этого тоже не знаем. И ветвления тоже зависят от архитектуры. Есть несколько водопроводов с предсказанием переходаи или нет?


yacc> Тогда самый быстрый способ - динамическая подлинковка.

Зобудь о ней. Это тоже вызов. Дорогой, с запоминанием точки возравта, со сбросом pipeline (привет ветвлениям).

yacc> самое быстрое вычисление.
Нет.
Самые быстрые — это линейные, когда нет переходов вообще, нет вызовов, линейное и водопровод разогнана, всё в кэше сидит. Это не считая математических преоброзований.
   31.031.0

AXT

инженер вольнодумец
★★
AXT> Сейчас уже не имеет смысла — ещё летом всё было заборото оптимизацией "в лоб", так что интересно только в смысле поиграться с математикой. Матрицы в стандарте ITU-T H.264, если нет доступа к, то могу выложить матрицы.

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

29557484
74740-74
84-29-7455
55-8474-29


Это самая частая, есть ещё:

64646464
8336-36-83
64-64-6464
36-8383-36


Ну тут структура очевидна ...

PS: Я монстр, я таки обманул тег csv :) Но только в предварительном просмотре :(
PPS: Таки есть режим auto :)
   13.0.782.22013.0.782.220
Это сообщение редактировалось 26.12.2014 в 17:52

yacc

старожил
★★★
Mishka> И последовательно лежит — это как?
Это с точки зрения адресации процессором.
По столбцам или по строкам - это уже компилятор преобразует [i][j] в соответствующий код

Mishka> И записывать обратно в вектор не можем до окончания операции полностью, т.к.о значение последнего элементо вектора-результата зависит от каждого элемента входного вектора.
В выходной? Это почему? Там операции вида out[k] += in[i] * matix[i][j] так и просятся.
Но быстрее использовать регистровую переменную.

Mishka> Зобудь о ней. Это тоже вызов. Дорогой, с запоминанием точки возравта, со сбросом pipeline (привет ветвлениям).
Ты меня не понял - я имею ввиду загрузку DLL с отображением ф-ии calc в память. А в ф-ии матрица самим набором вычислений представлена. Ее код генерирует отдельная программа-препроцессор.

Mishka> Самые быстрые — это линейные, когда нет переходов вообще, нет вызовов, линейное и водопровод разогнана, всё в кэше сидит. Это не считая математических преоброзований.
Там и нет переходов.
Там операции вида ( псевдо ассемблер, заглавными буквами - КОНСТАНТЫ )
- обнулить регистр С
- загрузить в регистр Б из памяти in + OFFSET1 ( например mov bx, _in + 25 )
- загрузить в регистр А константу MULT1
- умножить регистр А на регистр Б с результатом в А
- сложить регистр А с регистром С с результатов в А
- копировать А в С
...
- копировать С в область памяти out + OFFSET2

и НИКАКИХ ветвлений.

Надо несколько константных матриц? - загрузи несколько DLL вначале в ф-ии calc1, ... calcN
а потом дергай по-усмотрению.
Загрузка делается один раз.
Если вообще особо не меняется - можно и статически подлинковать эти ф-ии в исполняемый файл.
   
Это сообщение редактировалось 26.12.2014 в 14:49

AXT

инженер вольнодумец
★★
Mishka> Самые быстрые — это линейные, когда нет переходов вообще,

У некоторых DSP для этого имеется мощная фишка — аппаратные циклы. То есть, например (бывают варианты) цикл устроен так (псевдоассемблер):

R1 = (P0)+;
R2 = (P1)+;
LOOP (R0);
R3 = R1 + R2 || R1 = (P0)+;
MNOP || R2 = (P1)+ || (P2)+ = R3;
ENDLOOP;

В цикле автоматом выбираются и исполняются только строчки между LOOP и ENDLOOP, т.е. складываем значения за два такта. Счётчик цикла обрабатывается забесплатно.
   13.0.782.22013.0.782.220

Mishka

модератор
★★★
yacc> Это с точки зрения адресации процессором.
Копирование займёт время (и опять двойной вложенный цикл). Пробег по этим двум векторам — ещё цикл. И ты говоришь о том, что самое дорогое — это ветвления. Код, написанный ручками для каждого элемента — ни одного ветвления, столько же умножений и сложений с возможным уменьшеием слагаемых (если были нули в матрице), умножения (если были единицы и минус единицы), заменой умножения на сложение для небольших значений (если умножение медленное).

yacc> В выходной? Это почему? Там операции вида out[k] += in[i] * matix[i][j] так и просятся.
В исходный, как было одно время очень популярно на том FORTAN-е.

yacc> Ты меня не понял - я имею ввиду загрузку DLL с отображением ф-ии calc в память. А в ф-ии матрица самим набором вычислений представлена. Ее код генерирует отдельная программа-препроцессор.

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

yacc> Там и нет переходов.
Вызовы. Циклы.
   17.017.0

Mishka

модератор
★★★
AXT> В цикле автоматом выбираются и исполняются только строчки между LOOP и ENDLOOP, т.е. складываем значения за два такта. Счётчик цикла обрабатывается забесплатно.

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

в начало страницы | новое
 
Поиск
Настройки
Твиттер сайта
Статистика
Рейтинг@Mail.ru