Filter löschen
Filter löschen

Plotting wave solution at specific time

2 Ansichten (letzte 30 Tage)
bml727
bml727 am 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

Antworten (0)

Produkte


Version

R2019b

Community Treasure Hunt

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

Start Hunting!

Translated by