Clear Filters
Clear Filters

Plotting wave solution at specific time

1 view (last 30 days)
bml727
bml727 on 16 Apr 2020
I solved an equation u(x,t) and need to plot the solution at t=0.3, but I am not sure how to do it. I saved the time iteration and found that t=0.3 occurs at tplot(6). Any help would be great!
clear;
clc;
%% Problem 2
xstep = 0.1;
tstep = 0.05;
xstep2 = xstep*xstep;
tstep2 = tstep*tstep;
alpha = 2;
alpha2 = alpha*alpha;
lambda2 = alpha2*tstep2/xstep2;
xdomain = [0 1];
tdomain = [0 1];
nx = round((xdomain(2)-xdomain(1))/xstep);
nt = round((tdomain(2)-tdomain(1))/tstep);
xt0 = zeros((nx+1),1); % initial condition
dxdt0 = zeros((nx+1),1); % initial derivative
xold = zeros((nx+1),1); % solution at timestep k
x2old = zeros((nx+1),1); % solution at timestep k-1
xnew = zeros((nx+1),1); % solution at timestep k+1
% initial condition
pi = acos(-1.0);
for i=1:nx+1
xi = (i-1)*xstep;
if(xi>=0 && xi<=1)
xt0(i) = sin(2*pi*xi);
dxdt0(i) = alpha*pi*sin(2*pi*xi);
xold(i) = xt0(i)+dxdt0(i)*tstep;
xold(i) = xold(i) - 4*pi*pi*sin(2*pi*xi)*tstep2*alpha2;
end
end
x2old = xt0;
close all
x=linspace(xdomain(1),xdomain(2),nx+1);
%t=linspace(tdomain(1),tdomain(2),nx+1);
analy= sin(2*pi.*x).*(sin(4*pi.*0.3+cos(4*pi.*0.3)));
h1=plot(x,analy);
hold on;
h2=plot(x,xt0,'linewidth',2);
hold on;
h3=plot(x,xnew);
legend('Analytical','Initial','Final')
xlabel('x [m]');
ylabel('Displacement [m]');
set(gca,'FontSize',16);
tplot=zeros(1,nt);
for k=1:nt
time = k*tstep;
tplot(k)=time;
%tplot(6)=0.3
for i=1:nx+1
% Use periodic boundary condition, u(nx+1)=u(1)
if(i==1)
xnew(i) = 2*(1-lambda2)*xold(i) + lambda2*(xold(i+1)+xold(nx+1)) - x2old(i);
elseif(i==nx+1)
xnew(i) = 2*(1-lambda2)*xold(i) + lambda2*(xold(1)+xold(i-1)) - x2old(i);
else
xnew(i) = 2*(1-lambda2)*xold(i) + lambda2*(xold(i+1)+xold(i-1)) - x2old(i);
end
end
x2old=xold;
xold = xnew;
% if(mod(k,2)==0)
% h3.YData = xnew;
% refreshdata(h3);
% pause(0.5);
% end
end

Answers (0)

Categories

Find more on Graphics in Help Center and File Exchange

Products


Release

R2019b

Community Treasure Hunt

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

Start Hunting!