clear all; close all; %c=.75; c=1; delx=.1; a=-1; b=5; x=a:delx:b; J=(b-a)/delx; % set number of spatial steps. The nr. of dividing points is J+1 delt=.1; r=c*delt/delx; % c/rho. nsteps=20; u=f(x); % initialize the vector u. % Forward difference figure for n=1:nsteps v=u; for j=1:J u(j)=(1+r)*v(j)-r*v(j+1); end; if(n==1) plot(x,u,'g-+'); end; hold on; if(n==nsteps) plot(x,u,'r-s') hold on plot(x,sol(x,delt*nsteps)); % plot the solution with the method of characteristics %xlim([0 2]) %ylim([-1 1]) end; end; title('The solution using a forward difference approx.') xlabel('x','FontSize',14); ylabel('u_{fwd}(x,t)','FontSize',14); legend('n=1','n=nsteps','u_{exact}(x,4)') hold off; %break; u=f(x); % Upwind scheme figure for n=1:nsteps v=u; for j=2:J+1 u(j)=(1-r)*v(j)+r*v(j-1); end; if(n==1) plot(x,u,'g-+'); end; hold on; if(n==nsteps) plot(x,u,'r-s') hold on plot(x,sol(x,delt*nsteps)); % plot the solution with the method of characteristics end; end; title('The solution using an upwind difference approx.') xlabel('x','FontSize',14); ylabel('u_{uwd}(x,t)','FontSize',14); legend('n=1','n=nsteps','u_{exact}(x,4)') hold off; figure plot(x,f(x),'r'); % plot the initial solution hold on plot(x,sol(x,4)); % plot the solution with the method of characteristics hold off; title('Initial solution and Char.solution u(x(t),t)'); xlabel('x','FontSize',14); ylabel('u(x,t)','FontSize',14); legend('f','u_{char}(x,4)')