I need to plot of the time evolution of the wave equation in 1D when the Lax-Fri
ID: 3119041 • Letter: I
Question
I need to plot of the time evolution of the wave equation in 1D when the Lax-Friedrichs scheme, the Leapfrog scheme and the Lax-Wendroff scheme are used.I wrote MATLAB codes for these schemes but i can not find where i made the mistakes. I expect to plot same graphic that i've sent with the code .
%Gaussian
%The wave equation 4.2Lax-Fridrichs Scheme
clear all
close all
clc
xmin=0;
xmax=10;
N=100; %number of nodes-1
M=100; %number of nodes for time
tmin=0; %initial time
tmax=10; %max value of time
v=1.0 %input('velocity=') %velocity
%discretise the domain
dx=(xmax-xmin)/N;
x=xmin-dx:dx:xmax+dx; %ghost nodes
dt=(tmax-tmin)/M;
t=tmin-dt:dt:tmax+dt;
cfl=abs(v)*dt/dx
alfa=v*dt/dx;
sigma=1 %input('sigma=') %variance
mu=5 %input('mu=')
%initial gaussian
U(:,1)=(1/(sigma*(2*pi).^(1/2)))*exp(-((x-v*t-mu).^2)/(2*sigma.^2));
R(:,1)=-((v*(-mu-v*t+x).*exp(-((-mu-v*t+x).^2/(2*sigma.^2))))/(((2*pi).^(1/2))*(sigma.^3))); %v*du/dx %pde respect to x
S(:,1)=(v*(-mu-v*t+x).*exp(-((-mu-v*t+x).^2/(2*sigma.^2))))/(((2*pi).^(1/2))*(sigma.^3)); %du/dt % pde respect to t
for n=1:M
%calculate the scheme
for j=2:1:N+2
R(j,n+1)=0.5*(R(j+1,n)+R(j-1,n))+(alfa/2)*(S(j+1,n)-S(j-1,n));
S(j,n+1)=0.5*(S(j+1,n)+S(j-1,n))+(alfa/2)*(R(j+1,n)-R(j-1,n));
U(j,n+1)=U(j,n)+dt*S(j,n);
end
%update boundary condition
U(1,n+1)=U(N+1,n+1);
plot(x,U(:,n), '-o','markerfacecolor','b')
xlabel('x')
ylabel('U(x,t)')
%title(sprintf('time=%1.3f',t))
grid on
shg
pause(dt)
end
Explanation / Answer
I have resend the matlab code for the above mentioned wave equation .Please try with these codes.
%Gaussian
%The wave equation 4.2Lax-Fridrichs Scheme
clear all
close all
clc
xmin=0;
xmax=10;
N=100; %number of nodes-1
M=100; %number of nodes for time
tmin=0; %initial time
tmax=10; %max value of time
v=1.0 %input('velocity=') %velocity
%discretise the domain
dx=(xmax-xmin)/N;
x=xmin-dx:dx:xmax+dx; %ghost nodes
dt=(tmax-tmin)/M;
t=tmin-dt:dt:tmax+dt;
cfl=abs(v)*dt/dx
alfa=v*dt/dx;
sigma=1 %input('sigma=') %variance
mu=5 %input('mu=')
%initial gaussian
U(:,1)=(1/(sigma*(2*pi).^(1/2)))*exp(-((x-v*t-mu).^2)/(2*sigma.^2));
R(:,1)=-((v*(-mu-v*t+x).*exp(-((-mu-v*t+x).^2/(2*sigma.^2))))/(((2*pi).^(1/2))*(sigma.^3))); %v*du/dx %pde respect to x
S(:,1)=(v*(-mu-v*t+x).*exp(-((-mu-v*t+x).^2/(2*sigma.^2))))/(((2*pi).^(1/2))*(sigma.^3)); %du/dt % pde respect to t
for n=1:M
%calculate the scheme
for j=2:1:N+2
R(j,n+1)=0.5*(R(j+1,n)+R(j-1,n))+(alfa/2)*(S(j+1,n)-S(j-1,n));
S(j,n+1)=0.5*(S(j+1,n)+S(j-1,n))+(alfa/2)*(R(j+1,n)-R(j-1,n));
U(j,n+1)=U(j,n)+dt*S(j,n);
end
%update boundary condition
U(1,n+1)=U(N+1,n+1);
plot(x,U(:,n), '-o','markerfacecolor','b')
xlabel('x')
ylabel('U(x,t)')
%title(sprintf('time=%1.3f',t))
grid on
shg
pause(dt)
end