Plotting of error vs iteration number matlab

5 vues (au cours des 30 derniers jours)
Nima Vali
Nima Vali le 29 Oct 2020
Commenté : KSSV le 29 Oct 2020
I am going to plot the error as a function of iteration number in Jacoby calculation of second order diffrential equation. However, I am not getting any plot. Can you help. also, how I can get the plot at iteration number of 0.
clear all; clc ; close all;
%Setting up the parameters
dx=0.1 ; dy = dx ;
Lx=1 ; Ly =1;
x = (0:dx:Lx); y = (0:dy:Ly);
Nx=length(x); Ny=length(y) ;
if (Nx==Ny)
N=Nx;
end
[X Y] = meshgrid(x,y);
X=X' ; Y=Y' ;
%% calculating g(i,j)
g=-2*X.*(1-X)-2*Y.*(1-Y);
%% Analytical solution for tfinal
iterlimit=800;
f_exact=X.*(1-X).*Y.*(1-Y);
%%jacobi iterative method
f=zeros(Nx,Ny);
f_jacobi=zeros(Nx,Ny);
for n=1:iterlimit
for j=2:Ny-1
for i=2:Nx-1
f_jacobi(i,j)=(0.25)*(f(i,j+1)+f(i,j-1)+f(i+1,j)+f(i-1,j))-(0.25)*g(i,j)*dx^2 ;
end
end
f=f_jacobi;
B = reshape(f_jacobi-f_exact,[N*N,1]);
L1= norm(B,1)/numel(B);
plot(n,L1)
title ('L1 error vs iteration')
xlim([0,800])
end
  1 commentaire
KSSV
KSSV le 29 Oct 2020
You code has for loops which can be vectorised......think of your logic and code it again.

Connectez-vous pour commenter.

Réponse acceptée

KSSV
KSSV le 29 Oct 2020
clear all; clc ; close all;
%Setting up the parameters
dx=0.1 ; dy = dx ;
Lx=1 ; Ly =1;
x = (0:dx:Lx); y = (0:dy:Ly);
Nx=length(x); Ny=length(y) ;
if (Nx==Ny)
N=Nx;
end
[X Y] = meshgrid(x,y);
X=X' ; Y=Y' ;
%% calculating g(i,j)
g=-2*X.*(1-X)-2*Y.*(1-Y);
%% Analytical solution for tfinal
iterlimit=800;
f_exact=X.*(1-X).*Y.*(1-Y);
%%jacobi iterative method
f=zeros(Nx,Ny);
f_jacobi=zeros(Nx,Ny);
L1 = zeros(iterlimit) ;
for n=1:iterlimit
for j=2:Ny-1
for i=2:Nx-1
f_jacobi(i,j)=(0.25)*(f(i,j+1)+f(i,j-1)+f(i+1,j)+f(i-1,j))-(0.25)*g(i,j)*dx^2 ;
end
end
f=f_jacobi;
B = reshape(f_jacobi-f_exact,[N*N,1]);
L1(n)= norm(B,1)/numel(B);
end
plot(1:iterlimit,L1,'*r')
title ('L1 error vs iteration')
xlim([0,800])

Plus de réponses (0)

Catégories

En savoir plus sur Loops and Conditional Statements dans Help Center et File Exchange

Produits


Version

R2019b

Community Treasure Hunt

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

Start Hunting!

Translated by