clear;

n=5;
h=1/(n+1);

elem=((0:(n+1))')*h;

%% Scheme S1

maind=diag(ones(n+1,1));
upperd=diag(ones(n,1),-1);
lowerd=diag(ones(n,1),1);

A=(-lowerd+2*maind-upperd);
A(n+1,n+1)=1;

S1elem=((1:(n+1))')*h;
b=rhs214(S1elem)*h^2;
b(n+1)=h;

S1sol=A\b;
S1sol=[0; S1sol];

%% Scheme S2

mdiag=diag(ones(n+2,1));
updiag=diag(ones(n+1,1),1);
lwdiag=diag(ones(n+1,1),-1);

C=(-lwdiag+2*mdiag-updiag);
C(n+2,n+2)=1;
C(n+2,n+1)=0;
C(n+2,n)=-1;

S2elem=((1:(n+2))')*h;
d=rhs214(S2elem)*h^2;
d(n+2)=2*h;

S2sol=C\d;
S2sol=[0;S2sol(1:n+1)];

plot((0:n+1)'*h,S1sol,'g');
hold on;
plot((0:n+1)'*h,S2sol,'k');

N=15;
H=1/(N+1);
nodes=((0:(N+1))')*H;
solex=exsol214(nodes);
plot((0:N+1)'*H,solex,'ro');
hold off;


