Plotting of error vs iteration number matlab

3 ビュー (過去 30 日間)
Nima Vali
Nima Vali 2020 年 10 月 29 日
コメント済み: KSSV 2020 年 10 月 29 日
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 件のコメント
KSSV
KSSV 2020 年 10 月 29 日
You code has for loops which can be vectorised......think of your logic and code it again.

サインインしてコメントする。

採用された回答

KSSV
KSSV 2020 年 10 月 29 日
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])

その他の回答 (0 件)

カテゴリ

Help Center および File ExchangeLoops and Conditional Statements についてさらに検索

製品


リリース

R2019b

Community Treasure Hunt

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

Start Hunting!

Translated by