Trying to plot number of iterations vs gridsize for steady state

2 vues (au cours des 30 derniers jours)
Jayden Cavanagh
Jayden Cavanagh le 5 Juin 2021
I need to compute the number of iterations taken to reach the steady-state within a given tolerance. I then need to plot the final number of iterations against the grid size but I cannot for the life of me work out how I am supposed to do that. How am I suppose to set up and then plot those two things?
Code:
n=2;
nx=2^n;
nz=nx;
a=25;
b=25;
x = linspace(0, a, nx);
z = linspace(0, b, nz);
[X, Z] = meshgrid(x,z);
Tnp1 = zeros(nx, nz);
Tnp1(:,1) = 20;
Tnp1(:,end) = 20;
Tnp1(1,:) = 20+380*sin((X.*pi)/25)+205*sin((X.*5*pi)/25);
Tnp1(end,:) = 20;
err = 1;
tol = 1e-8;
k=0;
while err > tol
Tn = Tnp1;
k=k+1;
for i = 2:nx-1
for j = 2:nz-1
Tnp1(i,j) = (1/4)*(Tn(i+1,j)+Tn(i-1,j)+Tn(i,j+1)+Tn(i,j-1));
end
end
err = max(abs(Tnp1(:) - Tn(:)));
end
T = Tnp1;
plot(X,k)
end
end

Réponses (1)

SALAH ALRABEEI
SALAH ALRABEEI le 5 Juin 2021
Modifié(e) : SALAH ALRABEEI le 5 Juin 2021
Hi, this looks like a diffusion equation solved by FDM. you had to mistakes, I already corrected them, n is the step size or grid size, and the u should have used x not X in the boundary condition conditions. See the code below
clear N = [2,4]; for kk = 1: length(N) n=N(kk); nx=2^n; nz=nx; a=25; b=25; x = linspace(0, a, nx); z = linspace(0, b, nz); [X, Z] = meshgrid(x,z); Tnp1 = zeros(nx, nz); Tnp1(:,1) = 20; Tnp1(:,end) = 20; Tnp1(1,:) = 20+380*sin((x.*pi)/25)+205*sin((x.*5*pi)/25); Tnp1(end,:) = 20; err = 1; tol = 1e-8; k=0; while err > tol Tn = Tnp1; k=k+1; for i = 2:nx-1 for j = 2:nz-1 Tnp1(i,j) = (1/4)*(Tn(i+1,j)+Tn(i-1,j)+Tn(i,j+1)+Tn(i,j-1)); end end err = max(abs(Tnp1(:) - Tn(:))); end T = Tnp1; %subplot(2,2,kk) %surf(X,Z,T),grid %title([' Profile at grid = ' num2str(kk)])
Num_iter(kk) = k; end
figure(2) plot(N,Num_iter),grid
  2 commentaires
Jayden Cavanagh
Jayden Cavanagh le 5 Juin 2021
Modifié(e) : Jayden Cavanagh le 5 Juin 2021
It is a steady-state temperature field in a 2D rectangle solved using G-S method. So it isn't exactly what I am looking for. For my problem I have to choose the number of divisions in the domain to be ndiv = 2^n for n = 2 : 9 in both the x and z directions. n=2 was just a way to reduce processing power when testing. While your code works for N=[2,4] I get error code
Requested array exceeds the maximum possible variable size.
Error in linspace (line 44)
y = d1 + (0:n1).*(d2 - d1)./n1;
Error in hope (line 9)
x = linspace(0, a, nx);
when I try and run it for the maximum value of n=8 and nx=2^8. Any way the code could be changed to be able to composate for that?
SALAH ALRABEEI
SALAH ALRABEEI le 5 Juin 2021
Modifié(e) : SALAH ALRABEEI le 5 Juin 2021
%
clear
N = [2,4];
for kk = 1: length(N)
n=N(kk);
nx=2^n;
nz=nx;
a=25;
b=25;
x = linspace(0, a, nx);
z = linspace(0, b, nz);
[X, Z] = meshgrid(x,z);
Tnp1 = zeros(nx, nz);
Tnp1(:,1) = 20;
Tnp1(:,end) = 20;
Tnp1(1,:) = 20+380*sin((x.*pi)/25)+205*sin((x.*5*pi)/25);
Tnp1(end,:) = 20;
err = 1;
tol = 1e-8;
k=0;
while err > tol
Tn = Tnp1;
k=k+1;
for i = 2:nx-1
for j = 2:nz-1
Tnp1(i,j) = (1/4)*(Tn(i+1,j)+Tn(i-1,j)+Tn(i,j+1)+Tn(i,j-1));
end
end
err = max(abs(Tnp1(:) - Tn(:)));
end
T = Tnp1;
%subplot(2,2,kk)
%surf(X,Z,T),grid
%title([' Profile at grid = ' num2str(kk)])
Num_iter(kk) = k;
end
figure(2)
plot(N,Num_iter),grid

Connectez-vous pour commenter.

Catégories

En savoir plus sur Programming 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