% This is a script to test the theta method.

clear all;

theta=.5;
J=100;
T=2*pi;
N=1000;
dt=T/N;
dx=1/J;
b=0.01;
mu=b*dt/(dx*dx);
A1=mkA1(J,mu,theta);
A2=mkA2(J,mu,theta);
Yn=[1:J-1]'*0.0;
B=[1:J-1]'*0.0;
for n=0:N-1
    tn=n*dt;
    B(1)=mu*(1-theta)*sin(tn)+mu*theta*sin(tn+dt);
    Ynp1=A1\(A2*Yn+B);
    Yn=Ynp1;
end
U(J+1)=0;
U(1)=sin(T);
U(2:J)=Yn;
for j=0:J
   Xj(j+1)=dx*j;
end
