% zweite Version mit Matrizendarstellung
clear all
close all

I   = [0 1];
DIM = 40;

a0 = 1; a1 = 0; a2 = -1;                 % ODE: a2*f''+a1*f'+a0*f
rhs = @(x) (1+pi^2)*sin(pi*x);                  % right hand side
fl = 0; fr = 0;
%----------------------------------------------------------------
h = (I(2)-I(1))/(DIM-1);
x = I(1):h:I(2);
% rechte Seite
g = rhs(x);
% DGL in Matrix-Vektor-Darstellung
D0 = DiffMat(0,h,DIM);
D1 = DiffMat(1,h,DIM);
D2 = DiffMat(2,h,DIM);

A = a2*D2+a1*D1+a0*D0;
% Randwerte (Dirichlet) einbauen
A(1,1) = 1; g(1) = fl;
A(DIM,DIM) = 1; g(DIM) = fr;
% Loesung des LGS A*fh=g berechnen
fh = A\g';

subplot(1,2,1)
plot(x,fh,'k-');
grid on
subplot(1,2,2)
plot(x,abs(sin(pi*x)'-fh),'r-');
grid on
