Мне нужно повторно вычислить f(x)=exp(A*x) для крошечного переменного вектора-столбца x и огромной постоянной матрицы A (много строк, мало столбцов). Другими словами, x мало, а A*x много. Размеры моей проблемы таковы, что A*x занимает примерно столько же времени выполнения, сколько и часть exp().
Помимо расширения Тейлора и предварительного расчета диапазона значений exp(y) (при условии, что известен диапазон y значений A*x), который мне не удалось значительно ускорить (при сохранении точности) по сравнению с тем, что MATLAB делает сам по себе , Я думаю об аналитической переформулировке задачи, чтобы иметь возможность предварительно вычислить некоторые значения.
Например, я обнаружил, что exp(A*x)_i = exp(\sum_j A_ij x_j) = \prod_j exp(A_ij x_j) = \prod_j exp(A_ij)^x_j
Это позволило бы мне предварительно вычислить exp(A) один раз, но требуемое возведение в степень в цикле так же затратно, как и исходный вызов функции exp(), а умножения (\prod) должны выполняться дополнительно.
Есть ли какая-нибудь другая идея, которой я мог бы следовать, или решения в MATLAB, которые я мог пропустить?
Изменить: некоторые подробности
A имеет размер 26873856 на 81 (да, он такой огромный), поэтому x равен 81 на 1. nnz(A) / numel(A) равен 0.0012, nnz(A*x) / numel(A*x) равен 0.0075. Я уже использую разреженную матрицу для представления A, однако exp() разреженной матрицы больше не является разреженной. Так что на самом деле я сохраняю x неразреженным и вычисляю exp(full(A*x)), которое оказалось таким же быстрым/медленным, как full(exp(A*x)) (я думаю, что A*x в любом случае не является разреженным, так как x не является разреженным.) exp(full(A*sparse(x))) — это способ получить редкий A*x, но медленнее. Еще более медленные варианты — это exp(A*sparse(x)) (с удвоенной нагрузкой на память для неразреженной матрицы разреженного типа) и full(exp(A*sparse(x)) (что снова дает неразреженный результат).
sx = sparse(x);
tic, for i = 1 : 10, exp(full(A*x)); end, toc
tic, for i = 1 : 10, full(exp(A*x)); end, toc
tic, for i = 1 : 10, exp(full(A*sx)); end, toc
tic, for i = 1 : 10, exp(A*sx); end, toc
tic, for i = 1 : 10, full(exp(A*sx)); end, toc
Elapsed time is 1.485935 seconds.
Elapsed time is 1.511304 seconds.
Elapsed time is 2.060104 seconds.
Elapsed time is 3.194711 seconds.
Elapsed time is 4.534749 seconds.
Да, я рассчитываю опыт поэлементно, я обновляю приведенное выше уравнение, чтобы отразить это.
Еще одно редактирование: я пытался быть умным, но без особого успеха:
tic, for i = 1 : 10, B = exp(A*x); end, toc
tic, for i = 1 : 10, C = 1 + full(spfun(@(x) exp(x) - 1, A * sx)); end, toc
tic, for i = 1 : 10, D = 1 + full(spfun(@(x) exp(x) - 1, A * x)); end, toc
tic, for i = 1 : 10, E = 1 + full(spfun(@(x) exp(x) - 1, sparse(A * x))); end, toc
tic, for i = 1 : 10, F = 1 + spfun(@(x) exp(x) - 1, A * sx); end, toc
tic, for i = 1 : 10, G = 1 + spfun(@(x) exp(x) - 1, A * x); end, toc
tic, for i = 1 : 10, H = 1 + spfun(@(x) exp(x) - 1, sparse(A * x)); end, toc
Elapsed time is 1.490776 seconds.
Elapsed time is 2.031305 seconds.
Elapsed time is 2.743365 seconds.
Elapsed time is 2.818630 seconds.
Elapsed time is 2.176082 seconds.
Elapsed time is 2.779800 seconds.
Elapsed time is 2.900107 seconds.
A = rand(10000, 10). ИспользованиеexpA = exp(A); prod(bsxfun(@power, expA, x'), 2)действительно намного медленнее. Выполнениеexp(A*x)в функции mex тоже мало помогает, даже при использовании одинарной точности. Каков размерAдля вашего использования? - person buzjwa schedule 24.04.2014xявляется скаляром, аAвектором. Я не вижу никакого способа, как это можно ускорить. - person Bas Swinckels schedule 24.04.2014xменяется, аAнет? это происходит в шагах компрессионного зондирования? - person bla schedule 25.04.2014y=exp(A*x) .* x. Правильно или неправильно поняли? - person Horst Grünbusch schedule 29.04.2014