I am trying to solve the 2D Wave Equation numerically and plot the results.

14 次查看(过去 30 天)
Using the finite difference method, I have solved the 2D wave equation to get the following equation. When I try to plot the results, nothing happens.
The source function is working but plotting p(x,y,t) does not work. Here is the code:
%Numerically solve 2D Wave Equation
clear
clc
close all
%define constants
A = 2;
Omega = 2*pi/2.5;
StandardDev = 1;
speed = 1;
Nt = 100; %Number of temporal steps
Nx = 100; %Number of steps in x direction
Ny = 100; %Number of steps in y direction
tMin = 0;
tMax = 10;
xMin = 0;
xMax = 10;
xNot = xMax/2;
yMin = 0;
yMax = 10;
yNot = yMax/2;
t = linspace(tMin, tMax, Nt);
x = linspace(xMin, xMax, Nx);
y = linspace(yMin, yMax, Ny);
%Compute secondary parameters
DeltaT = t(2) - t(1);
DeltaX = x(2) - x(1);
DeltaY = y(2) - y(1);
%%Numerically solve PDE
%Initialize containers to hold all P values
f = zeros(Nx, Ny, Nt);
p = zeros(Nx, Ny, Nt);
%Set Initial Conditions
p(:,:,1) = 0;
p(:,:,2) = 0;
%Set Boundary Conditions
p(1,:,:) = 0;
p(Nx,:,:) = 0;
p(:,1,:) = 0;
p(:,Ny,:) = 0;
%Iterate through t values
for l=1:Nt
%Iterate through x values
for i=1:Nx
%Iterate through y values
for j=3:Ny
f(i,j,l) = A*sin(Omega*t(l))*exp(-1/(2*StandardDev^2)*((x(i) - xNot)^2 + (y(j) - yNot)^2));
end
end
end
for l=3:Nt
%Iterate through x values
for i=2:Nx
%Iterate through y values
for j=2:Ny
p(i,j,l) = 2*p(i,j,l-1) - p(i,j,l-2) + DeltaT^2*(f(i,j,l) + speed^2*((p(i+1,j,l-1) - 2*p(i,j,l-1) + p(i-1,j,l-1))/deltaX^2 + (p(i,j+1,l-1) - 2*p(i,j,l-1) + p(i,j-1,l-1))/deltaY^2));
end
end
end
%Plot the results
figure
for t=1:100
hold on
clf
surf(x, y, p(:,:,t))
zlim([-2 2])
colormap winter(1)
drawnow
end
hold off

回答(1 个)

Nikhil Sonavane
Nikhil Sonavane 2020-7-14
编辑:Nikhil Sonavane 2020-7-14
There are two spell errors in the following lines of code-
p(i,j,l) = 2*p(i,j,l-1) - p(i,j,l-2) + DeltaT^2*(f(i,j,l) + speed^2*((p(i+1,j,l-1) - 2*p(i,j,l-1) + p(i-1,j,l-1))/deltaX^2 + (p(i,j+1,l-1) - 2*p(i,j,l-1) + p(i,j-1,l-1))/deltaY^2));
deltaX should be DeltaX and similarly for deltaY as well. Even after correcting this there are some indexing errors in this line of code. In my suggestion you should correct this logic before moving forward.

类别

Help CenterFile Exchange 中查找有关 Geometry and Mesh 的更多信息

产品


版本

R2020a

Community Treasure Hunt

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

Start Hunting!

Translated by