Я кодирую решение уравнения Пуассона на двумерном прямоугольнике с использованием конечных элементов. Чтобы упростить код, я сохраняю дескрипторы базовых функций в массиве, а затем перебираю эти базовые функции, чтобы создать свою матрицу и правую часть. Проблема в том, что даже для очень грубых сеток это слишком медленно. Для сетки 9x9 (с использованием Dirichlet BC необходимо найти 49 узлов) это занимает около 20 секунд. Используя профиль, я заметил, что около половины времени тратится на доступ (а не на выполнение) моих базовых функций.
Профилировщик говорит matrix_assembly>@(x,y)bilinearBasisFunction(x,y,xc(k-1),xc(k),xc(k+1),yc(j-1),yc(j),yc(j+1)) (156800 calls, 11.558 sec)
, собственное время (без выполнения билинейного базового кода) составляет более 9 секунд. Есть идеи, почему это может быть так медленно?
Вот часть кода, при необходимости я могу опубликовать еще:
%% setting up the basis functions, storing them in cell array
basisFunctions = cell(nu, 1); %nu is #unknowns
i = 1;
for j = 2:length(yc) - 1
for k = 2:length(xc) - 1
basisFunctions{i} = @(x,y) bilinearBasisFunction(x,y, xc(k-1), xc(k),...
xc(k+1), yc(j-1), yc(j), yc(j+1)); %my code for bilinear basis functions
i = i+1;
end
end
%% Assemble matrices and RHS
M = zeros(nu,nu);
S = zeros(nu,nu);
F = zeros(nu, 1);
for iE = 1:ne
for iBF = 1:nu
[z1, dx1, dy1] = basisFunctions{iBF}(qx(iE), qy(iE));
F(iBF) = F(iBF) + z1*forcing_handle(qx(iE),qy(iE))/ae(iE);
for jBF = 1:nu
[z2, dx2, dy2] = basisFunctions{jBF}(qx(iE), qy(iE));
%M(iBF,jBF) = M(iBF,jBF) + z1*z2/ae(iE);
S(iBF,jBF) = S(iBF, jBF) + (dx1*dx2 + dy1*dy2)/ae(iE);
end
end
end
for
и изменяете значения параметров, используя каждый раз индексацию! И вы уверены, чтоnu
равно(length(xc)-2)*(length(yc)-2)
, чтобы избежать перераспределения в массивеbasisFunctions
ячеек? - person horchler   schedule 20.11.2013(jBF,qx(iE), qy(iE))
? - person Daniel   schedule 20.11.2013matrix_assembly
- это имя подфункции вbilinearBasisFunction
? - person horchler   schedule 20.11.2013matrix_assembly
и, может быть,bilinearBasisFunction
тогда. Даже если каждый вызов занимает несколько десятых миллисекунд, в сумме будет 156 800 вызовов функций. - person horchler   schedule 20.11.2013