clear all;
h = 0.01;
tau = 0.005;
x = 0:h:2;
t = 0:tau:1;
N = length(x);
M = length(t);
u = zeros([N,M]);
for m = 1:M
   u(1,m)= 0;
   u(N,m)= 0;
end
for n = 1:N
   u(n,1) = x(n)*sin(pi*x(n)^2);
   u(n,2) = u(n,1) + tau*x(n)^2*(x(n)-2)^2;
end
for m = 2:(M - 1)
   for n = 2:(N - 1)
       u(n,m+1) = 2*u(n,m) - u(n,m-1) + tau^2/h^2*(u(n+1,m)-2*u(n,m)+u(n-1,m));
   end
end
for m = 1:M
   plot(x,u(:,m))
   axis([0,2,-1,5])
   MM(m)= getframe;
end
movie(MM,3)