function u=fdtheta(u,x,h,T,lambda,theta); %%% Data: %%% u: initial data %%% x: spatial grid %%% h: cell size %%% T: final time %%% lambda, theta: characteristics of discretization %%% Finite Differences Theta Scheme for Parabolic Equation Nx=length(x)-1; dT=lambda*h^2; Nt=floor(T/dT); %%%% Construction of matrices A and B D=zeros(Nx-1,1); D(1)=1-2*(1-theta)*lambda; D(2)=(1-theta)*lambda; A=toeplitz(D); D(1)=1+2*theta*lambda; D(2)=-theta*lambda; B=toeplitz(D); for i=1:Nt u(2:Nx)=B\(A*u(2:Nx)); figure(1); plot(x,u); axis([-10 10 -1 1]); drawnow; end; [Nt*dT 1000*dT]