a = 1; b = 2;

xmany = linspace(a,b,500);
ymany = sqrt(xmany);
x = linspace(a,b,3);
y = sqrt(x);

lambda = baryzentrisch(x);
p = baryzentrischauswertung(x,xmany,y,lambda);
subplot(2,1,1);
plot(xmany,ymany,'b-',xmany,p,'r-');
makelarge
title(sprintf('Interpolation, Polynomgrad 3'));
subplot(2,1,2);
plot(xmany,(ymany-p),'r-');
makelarge
title('Interpolationsfehler');
disp('Max. Interpolationsfehler:');
max( abs(ymany-p) )

% Simpson-Regel
syms x;
exakt = 2/3* ( b^(3/2) - a^(3/2) );
quad  = (b-a) / 6 * ( y(1) + 4*y(2) + y(3) );
abs(quad - exakt)