Plotting of error vs iteration number matlab

3 Ansichten (letzte 30 Tage)
Nima Vali
Nima Vali am 29 Okt. 2020
Kommentiert: KSSV am 29 Okt. 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 Kommentar
KSSV
KSSV am 29 Okt. 2020
You code has for loops which can be vectorised......think of your logic and code it again.

Melden Sie sich an, um zu kommentieren.

Akzeptierte Antwort

KSSV
KSSV am 29 Okt. 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])

Weitere Antworten (0)

Kategorien

Mehr zu MATLAB finden Sie in Help Center und File Exchange

Produkte


Version

R2019b

Community Treasure Hunt

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

Start Hunting!

Translated by