How to convert to a parallel program?
Show older comments
n=41;
h=1/(n-1);
nt=40;
x(1)=0; y(1)=0;
for i=2:n
x(i)=x(i-1)+h;
y(i)=y(i-1)+h;
end
for j=1:n
for i=2:n
p0(i,j)=1;
end
end
tic
tau=0.01;
for jt=1:nt % Vaqt bo`yicha tsikl}
t(jt)=jt*tau;
for j=1:n
for i=1:n
pt(i,j)=exp(-2*x(i)*y(j)*t(jt));
end
end
% k+0.s cloy
for j=2:21 % Tenglamalar koeffitsientini hisoblash
for i=2:n-1
a(i)=1;
c(i)=1;
b(i)=2+2*h*h/tau;
ff=h*h*(-2*(2*t(jt)*t(jt)*(x(i)*x(i)+y(j)*y(j))+x(i)*y(j))*exp(-2*x(i)*y(j)*t(jt)));
d(i)=2*h*h/tau*p0(i,j)+(p0(i,j-1)-2*p0(i,j)+p0(i,j+1))+ff;
end
%Pragonka
aa1(1)=0; bb1(1)=1;
for i=2:n-1
aa1(i)=c(i)/(b(i)-aa1(i-1)*a(i));
bb1(i)=(d(i)+bb1(i-1)*a(i))/(b(i)-aa1(i-1)*a(i));
end
p(n,j)=exp(-2*y(j)*t(jt));
for i=n-1:-1:1
p(i,j)=aa1(i)*p(i+1,j)+bb1(i);
end
end %j
for j=40:-1:22 % Tenglamalar koeffitsientini hisoblash
for i=2:n-1
a(i)=1;
c(i)=1;
b(i)=2+2*h*h/tau;
ff=h*h*(-2*(2*t(jt)*t(jt)*(x(i)*x(i)+y(j)*y(j))+x(i)*y(j))*exp(-2*x(i)*y(j)*t(jt)));
d(i)=2*h*h/tau*p0(i,j)+(p0(i,j-1)-2*p0(i,j)+p0(i,j+1))+ff;
end
%Pragonka
aa1(1)=0; bb1(1)=1;
for i=2:n-1
aa1(i)=c(i)/(b(i)-aa1(i-1)*a(i));
bb1(i)=(d(i)+bb1(i-1)*a(i))/(b(i)-aa1(i-1)*a(i));
end
p(n,j)=exp(-2*y(j)*t(jt));
for i=n-1:-1:1
p(i,j)=aa1(i)*p(i+1,j)+bb1(i);
end
end %j
for i=1:n
p(i,n)=(4.0*p(i,n-1)-p(i,n-2))/3;
p(i,1)=(4.0*p(i,2)-p(i,3))/3;
end
for j=1:n
for i=1:n
p0(i,j)=p(i,j);
end
end
% k+1 sloy
for i=2:21 % Tenglamalar koeffitsientini hisoblash }
for j=2:n-1
a(j)=1;
c(j)=1;
b(j)=2+2*h*h/tau;
ff=h*h*(-2*(2*t(jt)*t(jt)*(x(i)*x(i)+y(j)*y(j))+x(i)*y(j))*exp(-2*x(i)*y(j)*t(jt)));
d(j)=2*h*h/tau*p0(i,j)+(p0(i-1,j)-2*p0(i,j)+p0(i+1,j))+ff;
end
%Pragonka
aa1(1)=0; bb1(1)=1;
for j=2:n-1
aa1(j)=c(j)/(b(j)-aa1(j-1)*a(j));
bb1(j)=(d(j)+bb1(j-1)*a(j))/(b(j)-aa1(j-1)*a(j));
end
p(i,n)=exp(-2*x(i)*t(jt));
for j=n-1:-1:1
p(i,j)=aa1(j)*p(i,j+1)+bb1(j);
end
end %i
for i=40:-1:22 % Tenglamalar koeffitsientini hisoblash }
for j=2:n-1
a(j)=1;
c(j)=1;
b(j)=2+2*h*h/tau;
ff=h*h*(-2*(2*t(jt)*t(jt)*(x(i)*x(i)+y(j)*y(j))+x(i)*y(j))*exp(-2*x(i)*y(j)*t(jt)));
d(j)=2*h*h/tau*p0(i,j)+(p0(i-1,j)-2*p0(i,j)+p0(i+1,j))+ff;
end
%Pragonka
aa1(1)=0; bb1(1)=1;
for j=2:n-1
aa1(j)=c(j)/(b(j)-aa1(j-1)*a(j));
bb1(j)=(d(j)+bb1(j-1)*a(j))/(b(j)-aa1(j-1)*a(j));
end
p(i,n)=exp(-2*x(i)*t(jt));
for j=n-1:-1:1
p(i,j)=aa1(j)*p(i,j+1)+bb1(j);
end
end %i
for j=1:n
p(n,j)=(4.0*p(n-1,j)-p(n-2,j))/3;
p(1,j)=(4.0*p(2,j)-p(3,j))/3;
end
for i=1:n
for j=1:n
p0(i,j)=p(i,j);
end
end
end % jt
toc
for i=1:n
for j=1:n
b(j)=(p(i,j)-pt(i,j));
max=b(j);
end
end
for j=1:n
px(j)=p(21,j);
pxt(j)=pt(21,j);
end
t
p
pt
b
max
4 Comments
Jan
on 17 Sep 2021
What does "parallel program" mean? Do you want to vectrorize the code to increase the speed?
We cannot run the code due to the missing variables nt, n, x, y.
The outer loop is useless, because the contents of pt is overwritten in each iteration. Then it is much cheaper to omit the outer loop and run the code for jt=nt only.
Do you pre-allocate pt before the loops? Iteratively growing arrays need an exponentially growing amount of RAM, even if the final output is not such large.
The current version of your code does not need loops:
pt = exp(-2 * x(:) .* y(:).' * (nt * tau));
% ^^^ ^^^ maybe this can be omitted depending on
% what x and y are
Shixnazar Ismailov
on 17 Sep 2021
Shixnazar Ismailov
on 17 Sep 2021
Accepted Answer
Jan
on 17 Sep 2021
Pre-allocation is essential. Do not let array grow iteratively, but create them with the final size:
n = 41;
h = 1 / (n-1);
nt = 40;
% x(1) = 0;
% y(1) = 0;
% Better:
x = zeros(1, n); % Pre-allocation
y = zeros(1, n);
for i = 2:n
x(i) = x(i - 1) + h;
y(i) = y(i - 1) + h;
end
In this case you do not even need a loop:
x = (0:n - 1) * h;
y = x;
Replace:
for j = 1:n
for i = 2:n
p0(i,j)=1;
end
end
by
p0 = ones(n, n);
This is faster and nicer. Pre-allocate t also:
tau = 0.01;
t = (1:nt) * tau;
for jt = 1:nt % Vaqt bo`yicha tsikl}
% pt without loops:
pt = exp(-2 * t(jt) * x.' .* y);
% Moved out of the loops:
a = ones(1, n-1);
c = ones(1, n-1);
b = repmat(2 + 2 * h * h / tau, 1, n-1);
% k+0.s cloy
for j = 2:21 % Tenglamalar koeffitsientini hisoblash
for i = 2:n-1
% Don't do this repeatedly inside the loops - see above:
% a(i) = 1;
% c(i) = 1;
% b(i) = 2 + 2 * h * h / tau;
% You have evaluated exp(-2*x(i)*y(j)*t(jt)) already above.
% Don't do it again.
ff = h * h * (-2 * (2 * t(jt)*t(jt) * (x(i)*x(i) + y(j)*y(j)) + x(i)*y(j)) * pt(i, j));
d(i) = 2 * h * h / tau * p0(i,j) + (p0(i,j-1) - 2*p0(i,j) + p0(i,j+1)) + ff;
end
%Pragonka
aa1 = zeros(1, n-1); % Pre-allocate!
bb1 = zeros(1, n-1); % Pre-allocate!
bb1(1) = 1;
for i = 2:n-1
aa1(i) = c(i) / (b(i) - aa1(i-1) * a(i));
bb1(i) = (d(i) + bb1(i-1) * a(i)) / (b(i) - aa1(i-1) * a(i));
end
p(n,j) = exp(-2 * y(j) * t(jt));
for i = n-1:-1:1
p(i,j) = aa1(i) * p(i+1,j) + bb1(i);
end
end
Replace:
for j=1:n
for i=1:n
p0(i,j) = p(i,j);
end
end
by
p0 = p;
Please try to implement the shown methods in your code.
Use spaces around operators and a proper indentation.
for i = 1:n
for j = 1:n
b(j) = p(i,j) - pt(i,j);
max = b(j); % ???
end
end
What is the purpose of this code? b(j) is overwritten in each iteration of i. Finally the value of max is p(n,n) - pt(n,n) and does not depend on the loop in any way. It is not worth to improve the speed of a code, when it does not compute anything useful.
1 Comment
Shikhnazar Ismailov
on 24 Feb 2022
function []=diagonal()
n=41; ns=(n+1)/2; h=1/(n-1); nt=360;
x = zeros(1, n);
y = zeros(1, n);
x = (0:n - 1) * h;
y = x;
p0 = ones(n);
tau=0.01;
for jt=1:nt
t(jt)=jt*tau;
pt = exp(-2 * t(jt) * x.' .* y);
end
tic
for jt=1:nt % Vaqt bo`yicha tsikl}
a = ones(1, n-1);
c = ones(1, n-1);
b = repmat(2 + 2 * h * h / tau, 1, n-1);
aa1 = zeros(1, n-1);
bb1 = zeros(1, n-1);
bb1(1) = 1;
for j=2:ns
parfor i=2:n-1
ff=h*h*(-2*(2*t(jt)*t(jt)*(x(i)*x(i)+y(j)*y(j))+x(i)*y(j))*exp(-2*x(i)*y(j)*t(jt)));
d=2*h*h/tau*p0(i,j)+(p0(i,j-1)-2*p0(i,j)+p0(i,j+1))+ff;
aa1(i)=c/(b-aa1(i-1)*a);
bb1(i)=(d+bb1(i-1)*a)/(b-aa1(i-1)*a);
end
p(n,j)=exp(-2*y(j)*t(jt));
for i=n-1:-1:1
p(i,j)=aa1(i)*p(i+1,j)+bb1(i);
end
end %j
for j=n-1:ns-1
parfor i=2:n-1
ff=h*h*(-2*(2*t(jt)*t(jt)*(x(i)*x(i)+y(j)*y(j))+x(i)*y(j))*exp(-2*x(i)*y(j)*t(jt)));
d=2*h*h/tau*p0(i,j)+(p0(i,j-1)-2*p0(i,j)+p0(i,j+1))+ff;
aa1(i)=c/(b-aa1(i-1)*a);
bb1(i)=(d+bb1(i-1)*a)/(b-aa1(i-1)*a);
end
p(n,j)=exp(-2*y(j)*t(jt));
for i=n-1:-1:1
p(i,j)=aa1(i)*p(i+1,j)+bb1(i);
end
end %j
for i=1:n
p(i,n)=(4.0*p(i,n-1)-p(i,n-2))/3;
p(i,1)=(4.0*p(i,2)-p(i,3))/3;
end
p0=p;
for i=2:ns
parfor j=2:n-1
ff=h*h*(-2*(2*t(jt)*t(jt)*(x(i)*x(i)+y(j)*y(j))+x(i)*y(j))*exp(-2*x(i)*y(j)*t(jt)));
d=2*h*h/tau*p0(i,j)+(p0(i-1,j)-2*p0(i,j)+p0(i+1,j))+ff;
aa1(j)=c/(b-aa1(j-1)*a);
bb1(j)=(d+bb1(j-1)*a)/(b-aa1(j-1)*a);
end
p(i,n)=exp(-2*x(i)*t(jt));
for j=n-1:-1:1
p(i,j)=aa1(j)*p(i,j+1)+bb1(j);
end
end %i
for i=n-1:ns-1
parfor j=2:n-1
ff=h*h*(-2*(2*t(jt)*t(jt)*(x(i)*x(i)+y(j)*y(j))+x(i)*y(j))*exp(-2*x(i)*y(j)*t(jt)));
d=2*h*h/tau*p0(i,j)+(p0(i-1,j)-2*p0(i,j)+p0(i+1,j))+ff;
aa1(j)=c/(b-aa1(j-1)*a);
bb1(j)=(d+bb1(j-1)*a)/(b-aa1(j-1)*a);
end
p(i,n)=exp(-2*x(i)*t(jt));
for j=n-1:-1:1
p(i,j)=aa1(j)*p(i,j+1)+bb1(j);
end
end %i
for j=1:n
p(n,j)=(4.0*p(n-1,j)-p(n-2,j))/3;
p(1,j)=(4.0*p(2,j)-p(3,j))/3;
end
p0=p;
end % jt
toc
px = p0(21,:);
pxt = pt(21,:);
x=0:h:1;
y=0:h:1;
figure ('Position',[600 60 700 600]);
subplot(3,2,1);
meshc(x,y,pt) %Uch o'lchovli 3D grafika
title('АНИК ЕЧИМ 3D ГРАФИГИ')
subplot(3,2,2);
meshc(x,y,p) %Uch o'lchovli 3D grafika
title('CОНЛИ ЕЧИМ 3D ГРАФИГИ')
subplot(3,2,3);
contour(x,y,pt,'ShowText','on','LineWidth',2); %kontur chizish
title('АНИК ЕЧИМ КОНТУР ГРАФИГИ')
subplot(3,2,4);
contour(x,y,p,'ShowText','on','LineWidth',2); %kontur chizish
title('CОНЛИ ЕЧИМ КОНТУР ГРАФИГИ')
subplot(3,2,5);
plot(x,px,x,pxt); %Ikki o'lchovli grafika
title('КЕСИМДА АНИК ВА CОНЛИ ЕЧИМ ЎЗГАРИШИ X БЎЙИЧА')
end
The variable 'aa1' in a parfor cannot be classified.
More Answers (0)
Categories
Find more on Loops and Conditional Statements in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
