How to plot 2D line graph to compare the approximate solution with the actual solution?

2 views (last 30 days)
I am trying to solve the following problem.
I wrote the following code to solve this problem. I stored my approximate solution in the matrix U of size Nx+2 times Ny+2. I want to compare my approximated solution with the actual solution. I am not sure how to plot a matrix in to 2D line graph. Can anyone please help me with this step?
%------Construction of meshgrids-------------------
Lx = 1;
Ly = 1;
%Ny=Nx;
Nx = 3; % Number of interior grid points in x direction
Ny = 3; % Number of interior grid points in y direction
dx = Lx/(Nx+1); % Spacing in the x direction
dy = Ly/(Ny+1); % Spacing in the y direction
x=0:dx:Lx;
y=0:dy:Ly;
[X,Y] = meshgrid(x,y);%2d array of x,y values
X = X'; % X(i,j), Y(i,j) are coordinates of (i,j) point
Y = Y';
% -------------------------------------------------------
Iint = 1:Nx+1; % Indices of interior point in x direction
Jint = 1:Ny+2; % Indices of interior point in y direction
Xint = X(Iint,Jint);
Yint = Y(Iint,Jint);
U = zeros(Nx+2,Ny+2); % Define U to store the answer
%--------------------------------------------------------
uinit = zeros(Nx+1,Ny+2);
u0 = @(x,y) cos(2*pi*x).*cos(2*pi*y);
U(1,:) = u0(x(1),y(Jint));
U(Nx+2,:)=U(1,:);
uinit(1,:)=U(1,:);
uinit(Nx,:)=uinit(1,:);
F1 = reshape(uinit, (Nx+1)*(Ny+2),1);
%---------------------------------------------
%assembly of the tridiagonal matrix(LHS)
sx = 1/(dx^2);
sy = 1/(dy^2);
e=ones(Nx,1);
A = zeros(Nx,Nx+2);
B = zeros(Nx+1,Nx+2);
T=spdiags([sx.*e,(-2*sx)+(-2*sy).*e,sx.*e],-1:1,Nx,Nx);
A(:,2:Nx+1) = T;
B(1:Nx,:)=A;
B(Nx+1,1)=-1;
B(Nx+1,2)=1;
B(Nx+1,Nx+1)=1;
B(Nx+1,Nx+2)=-1;
%T(1,Nx+1)= sx;
%T(Nx+1,1)= sx;
D = zeros(Nx+1,Nx+2);
D1 = zeros(Nx,Nx+2);
D2=spdiags(sy.*e,0,Nx,Nx);
D1(:,2:Nx+1)=D2;
D(1:Nx,:)=D1;
C=blktridiag(B,D,D,Ny+2);
for i=1:Nx
for j=1:Nx+1
C(i,j)=(1/2).*C(i,j);
C((Nx+1)*(Ny+1)+i,(Nx+2)*(Ny+1)+j)=(1/2).*C((Nx+1)*(Ny+1)+i,(Nx+2)*(Ny+1)+j);
end
end
%---------------------------------------------------------------
%----------Compuating the R.H.s term----------------------------
f = @(x,y) (-8*pi*pi).*cos(2*pi*x).*cos(2*pi*y);
%Solve the linear system
rhs = f(Xint,Yint);
rhs_new = zeros(Nx+1,Ny+2);
rhs_new(1:Nx,:)=rhs(2:Nx+1,:);
%for j=1:Ny+2
% for i=2:Nx+3
% rhs(i,j) = (-8*pi*pi)*cos(2*pi*x(i)).*cos(2*pi*y(j));
% end
%end
for i=1:Nx
rhs_new(i,1)=(1/2).*rhs_new(i,1);
rhs_new(i,Ny+2)=(1/2).*rhs_new(i,Ny+2);
end
%convert the rhs into column vector
F = reshape(rhs_new, (Nx+1)*(Ny+2),1);
F2 = F - sx.* F1;
uvec = C\F2;
v= reshape(uvec, Nx+2,Ny+2);
U(2:Nx+1,:)=v(2:Nx+1,:);
%----------Computation of Exact solution--------------------
ue = zeros(Nx+2,Ny+2);
ue(:,:) = cos(2*pi*X) .* cos(2*pi*Y); % Exact Solution
%Error
err = norm(U(:,:) - ue(:,:))/norm(ue(:,:));

Answers (1)

Divyam
Divyam on 5 Sep 2024
Hi @kaps,
You can plot a surface plot using the "surf" function to compare the approximate solution with the exact solution.
% Plot both solutions as surface plots
figure;
% Create a subplot with 1 row and 2 columns
subplot(1, 2, 1);
% Plot the approximate solution
surf(X, Y, U);
colorbar;
title('Approximate Solution (U)');
xlabel('x');
ylabel('y');
zlabel('U');
axis tight;
subplot(1, 2, 2);
% Plot the exact solution
surf(X, Y, ue);
colorbar;
title('Exact Solution (ue)');
xlabel('x');
ylabel('y');
zlabel('ue');
axis tight;
For more information regarding other methods to plot surface and mesh plots, refer to this documentation link: https://www.mathworks.com/help/matlab/surface-and-mesh-plots-1.html

Categories

Find more on Function Creation 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!