% TESTISTARAL  Test of Istaral
%
global alpha beta abj abj1 abj2 abj3 abj4 abl1 abl2 abl3 abl4
f0='%3.0f %3.0f %3.0f %15.12f %15.12f %10.2e %9.2e\n';
f1='%3.0f %3.0f %3.0f %15.12f %15.12f %10.2e\n';
disp(' al  be   n       Ialg            Ilog          erra      errl')
n=8;
for alpha=[0 3 5]
  for beta=[0 4 6]
    abj=r_jacobi01(n,alpha,beta); abj1=r_jacobi01(n,alpha,0);
    abj2=r_jacobi01(n,beta,0); abj3=r_jacobi01(n,alpha,2*beta+1);
    abj4=r_jacobi01(n,beta,2*alpha+1);
    abl1=r_jacobilog(n,alpha,0); abl2=r_jacobilog(n,beta,0);
    abl3=r_jacobilog(n,alpha,2*beta+1);
    abl4=r_jacobilog(n,beta,2*alpha+1);
    Istaralgex=(gamma(alpha+1)*gamma(beta+1)/gamma(alpha+beta+2))^2;
    Istarlogex=-3/2;
    for n=1:8
      [Istaralg,Istarlog]=Istaral(n);
      errIstaralg=abs((Istaralg-Istaralgex)/Istaralgex);
      if alpha==0 & beta==0
        errIstarlog=abs((Istarlog-Istarlogex)/Istarlogex);
        fprintf(f0,alpha,beta,n,Istaralg,Istarlog, ...
      errIstaralg,errIstarlog)
      else
        fprintf(f1,alpha,beta,n,Istaralg,Istarlog, ...
          errIstaralg)
      end
    end
    fprintf('\n')
  end
end
