Improve code for the numerical aproximation of the 1D fractional heat equation
4 次查看(过去 30 天)
显示 更早的评论
I need a way to improve the time it takes to calculate the matrix. Right now it takes around 19 seconds for the small interval [-100,100] and h = .125. and I need it to work for h = e-6 and and the interval [a,b] = [-5000,5000]
tic
%Bound
a = -100;
b = 100;
T = 1;
%Alpha value
alpha = 1;
h = 0.125;
tau = h^alpha/100;
%Boundary conditions
f1 = @(x) 1./(1+x.^2);
%Exact solution
y = @(x,t) (t+1)./((t+1).^2 + x.^2);
%Draw of the aproximation and the error
DFNLE(alpha,a,b,T,h,tau,f1,y)
With the main function being
function DFNLE(alpha,a,b,T,h,tau,f1,y)
n = (b-a)/h - 1;
n = floor(n);
x = a:h:b;
m = T / tau;
m = nearest(m);
t = 0:tau:T;
w = zeros(n+2,1);
v1 = zeros(n,1);
v2 = zeros(n,1);
%Restrain
theta = h^(1/3);
%Constant C
C = 2^alpha * gamma(1/2+alpha/2) / (sqrt(pi)*abs(gamma(-alpha/2)));
%Known values of u
u = zeros(n+2,m+1);
u(:,1) = f1(x);
%Loop to create the weights
for k = 1:n+2
if k*h > theta
w(k) = C*h^(-alpha)/k^(alpha+1);
end
end
%CFL condition
if tau > 1/sum(w)
fprintf('Not convergent tau = %s', tau);
return
end
% Main loop
for j = 1:m
for i = 1:n+2
for k = 1:n+2
if i+k <= n+2
v1(k) = (u(i+k,j)-u(i,j))*w(k);
else
v1(k) = -u(i,j)*w(k);
end
if i-k >= 1
v2(k) = (u(i-k,j)-u(i,j))*w(k);
else
v2(k) = -u(i,j)*w(k);
end
end
u(i,j+1) = u(i,j)+tau*(sum(v1)+sum(v2));
end
end
final = T;
[X,T] = meshgrid(x,t);
u = sparse(u);
%Error
err = max(max(abs(u' - y(X,T))));
fprintf('Error, epsilon = %d', err)
% Drawing
% surf(X,T,u');
% axis([-10 10 0 final 0 1])
% shading interp
% title(['Explicit method for alpha = ', num2str(alpha)])
% xlabel('x')
% ylabel('t')
% zlabel('function u')
end
0 个评论
回答(1 个)
Shadaab Siddiqie
2021-4-16
From my understanding you want to improve the performance of above code. MATLAB is optimized to work on matrices and vectorizing the code improves MATLAB performance compared to using nested FOR loops.
The following two examples use the functions ARRAYFUN and CELLFUN to compute the same final matrix without using FOR loops. When the size of the matrix is small, all these approaches are similar in terms of the time needed to execute them. However, when the matrix size increases, these approaches are much faster.
%% 1st approach: using arrays and ARRAYFUN
% Step 1: Compute the cumulative variance vector (last column of the target matrix Y).
tic
f = @(k) var(Z(1:k));
w = arrayfun(f, (1:numel(Z)).');
% Step 2: Assemble the matrix (which is built from zeros
% and subvectors of the cumulative variance vector).
Y3 = triu( w*ones(1, numel(w)) );
toc
%% 2nd approach: using cell arrays and CELLFUN
A = Z';
n = length(A);
tic
Y2=triu(repmat(cellfun(@(x) var(x(isfinite(x))), num2cell(repmat(A,n,1) + triu(NaN(n, n),1),2)),1,n));
toc
%% Comparison
isequal(Y, Y2, Y3)
You can incorporate these functions and vectorize your code to improve it's performance.
0 个评论
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Mathematics and Optimization 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!