%{ # Diffusion equation without time marching (but reversed) Here is like in [advection_global_reversed.m](), except that we do this for the diffusion equation (which should be much more ill-behaved, and more sensitive to errors). %} clear all; clf % parameters Nx=61; % number of gridpoints Nt=60; % gridpoints in time Lx=10; % domain length in x Lt=4; % duration in time xpos=5; % position of the final condition mu=0.1; % differentiation [d.x,d.xx,d.wx,x]=dif1D('cheb',0,Lx,Nx); [d.y,d.yy,d.wy,y]=dif1D('fd',0,Lt,Nt,5); [D,l,X,Y,Z,I,NN]=dif2D(d,x,y); % rename y to time derivative D.t=D.y; D.tt=D.yy; t=y; T=Y; l.left=[l.cbl;l.left;l.ctl]; l.right=[l.cbr;l.right;l.ctr]; l.start=l.bot; l.end=l.top; % system matrices A=D.t-mu*D.xx; b=zeros(NN,1); %{ Here this is the "final guess". %} % final guess f0=exp(-((X-xpos)/1).^2); %{ # Boundary conditions The "left", "right" and "end" boundaries should be like in the initial guess. %} % boundary conditions loc=[l.end; l.right; l.left]; C=I(loc,:); b(loc)=C*f0(:); A(loc,:)=C; %{ # Solve system %} f=A\b; f=reshape(f,Nt,Nx); % show evolution of f subplot(1,2,1); surf(X,T,f); view(2) xlabel('x'); ylabel('t'); title('time evolution of the string'); %{ # March in time of the initial condition Here this is some kind of validation, we use a different method to go back from the initial condition to the final one, and compare to the global solution. %} % march in time the initial condition I=eye(Nx); E=I; A=mu*d.xx; % boundary conditions loc=[1;Nx]; E(loc,:)=0; A(loc,:)=I(loc,:); % march in time matrix dt=t(2)-t(1); Mm=(E-A*dt/2); M=Mm\(E+A*dt/2); % initial condition q=f(1,:)'; % marching loop for ind=1:(Nt-1); q=M*q; end % plotting subplot(1,2,2); plot(x,q,'b',X(l.end),f(l.end),'r.-'); axis([0,Lx,0,1.2]) legend('march in time','global'); title('diffusion reversed') xlabel('x'); ylabel('final time'); set(gcf,'paperpositionmode','auto'); print('-dpng','-r80','diffusion_global_reversed.png'); %{ The subplot on the left shows the evolution in time of the solution. You see that this is ill-behaved because the solution becomes very large and oscillatory. ![Validation](diffusion_global_reversed.png) # Exercises/contributions * Please add exercices and contributions %}