%code1 for lecture 2; need u_0(x) and f(x, t)

    dt=input('  dt=  ');
    T=input(  'T=  ');

  for N=2:12
   
    d1=DM1(N);  d2=d1*d1;
    A=d2(2:N,2:N);

%initial condition 
    j=[1:1:N-1];  x=[cos(j*pi/N)]';  u=sin(pi*x);

     time=0.0;

     while time <= T 

      K1=A*u;
      K2=A*(u+0.5*dt*K1);
      K3=A*(u+0.5*dt*K2);
      K4=A*(u+dt*K3);
      u=u+dt*(K1+2*K2+2*K3+K4)/6;

      time=time+dt;
      
     end
%exact solution

       exact=exp(-pi^2*time)*sin(pi*x); error=abs(exact-u);

       xx(N-1)=N;
       err(N-1)=max(error);
       
       ss=max(error)

  end
  
  fprintf(1,      '%9.0f          %10.2e  \n', [xx;  err]);
