Трябва да изчисля 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)) (с удвоено въздействие върху паметта за неразредена матрица от тип sparse) и 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.
Да, изчислявам exp по елементи, актуализирам горното уравнение, за да отразя това.
Още една редакция: Опитах се да бъда умен, но с малък успех:
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