Improve code for the numerical aproximation of the 1D fractional heat equation

4 vues (au cours des 30 derniers jours)
SerchM
SerchM le 13 Avr 2021
Modifié(e) : SerchM le 16 Avr 2021
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

Réponses (1)

Shadaab Siddiqie
Shadaab Siddiqie le 16 Avr 2021
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.

Catégories

En savoir plus sur Sparse Matrices dans Help Center et File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by