% Script: TestCI2D
%
% Tests CubicInterp2D on two examples:
%      f=x^3 y - x y^3
%
%      f = sin(2 pi x) exp(-y^2)
%
n=20;
m=35;
a=-1;
b=1;
c=0;
d=3;
x=linspace(a,b,n);
y=linspace(c,d,m);
fA1=zeros(n,m);
fA2=zeros(n,m);
for i=1:n
   for j=1:m
      fA1(i,j)=y(j)*(x(i))^3 - x(i)*(y(j))^3;
      fA2(i,j)=sin(2*pi*x(i))*exp(-y(j)*y(j));
   end
end
xc1=0.273;
yc1=.9871;
xc2=-.9982;
yc2=2.981;
disp(' x^3 y - x y^3 : The interpolant should be exact ')
z=CubicInterp2D(xc1,yc1,a,b,c,d,fA1);
true=(xc1^3)*yc1-xc1*(yc1^3);
err1=abs(z-true);
z=CubicInterp2D(xc2,yc2,a,b,c,d,fA1);
err2=abs(z-((xc2^3)*yc2-xc2*(yc2^3)));
disp(sprintf('errors for middle and corner points %1.10f %1.10f',err1,err2))
%
disp(' sin(2*pi*x)*exp(-y^2) ')
z=CubicInterp2D(xc1,yc1,a,b,c,d,fA2);
err1=abs(z-(sin(2*pi*xc1)*exp(-yc1*yc1)));
z=CubicInterp2D(xc2,yc2,a,b,c,d,fA2);
err2=abs(z-(sin(2*pi*xc2)*exp(-yc2*yc2)));
disp(sprintf('errors for middle and corner points %1.10f %1.10f',err1,err2))
 