[lun,kprint,ipass]=pchqk5(lun,kprint,ipass);

function [lun,kprint,ipass]=pchqk5(lun,kprint,ipass);
persistent ansmlv bcoef d dcalc dermax derr err f fail fcalc fermax ferr firstCall i ierr ifail inbv j k knotyp n ndim nknots t termax terr tol tolz tsave work x zero ; if isempty(firstCall),firstCall=1;end; 

if isempty(i), i=0; end;
if isempty(ierr), ierr=0; end;
if isempty(ifail), ifail=0; end;
if isempty(inbv), inbv=0; end;
if isempty(j), j=0; end;
if isempty(knotyp), knotyp=0; end;
if isempty(k), k=0; end;
if isempty(ndim), ndim=0; end;
if isempty(nknots), nknots=0; end;
if isempty(n), n=9 ; end;
if isempty(bcoef), bcoef=zeros(1,2.*n); end;
if isempty(d), d=zeros(1,n); end;
if isempty(dcalc), dcalc=0; end;
if isempty(derr), derr=0; end;
if isempty(dermax), dermax=0; end;
if isempty(f), f=zeros(1,n); end;
if isempty(fcalc), fcalc=0; end;
if isempty(ferr), ferr=0; end;
if isempty(fermax), fermax=0; end;
if isempty(t), t=zeros(1,2.*n+4); end;
if isempty(terr), terr=0; end;
if isempty(termax), termax=0; end;
if isempty(tol), tol=0; end;
if isempty(tolz), tolz=0; end;
if isempty(tsave), tsave=zeros(1,2.*n+4); end;
if isempty(work), work=zeros(1,16.*n); end;
if isempty(x), x=zeros(1,n); end;
if isempty(zero), zero=0.0e0 ; end;
if isempty(fail), fail=false; end;
if isempty(ansmlv), ansmlv=0; end;
if isempty(err), err=0; end;
% relerr= @(err,ansmlv)  abs(err)/max(1.0e-5,abs(ansmlv));real :: relerr;
if firstCall,   x=[-2.2e0,-1.2e0,-1.0e0,-0.5e0,-0.01e0,0.5e0,1.0e0,2.0e0,2.2e0];  end;
if firstCall,   f=[0.0079e0,0.2369e0,0.3679e0,0.7788e0,0.9999e0,0.7788e0,0.3679e0,0.1083e0,0.0079e0];  end;
if firstCall,   d=[0.0000e0,0.3800e0,0.7173e0,0.5820e0,0.0177e0,-0.5696e0,-0.5135e0,-0.0778e0,-0.0025e0];  end;
firstCall=0;
relerr= @(err,ansmlv)  abs(err)./max(1.0e-5,abs(ansmlv));
ifail = 0;
tol = 100.*r1mach(4);
tolz = zero;
if( kprint>=3 )
writef(lun,['1', '\n ' , '\n ' ,repmat(' ',1,10),'TEST PCH TO B-SPLINE CONVERTER' ' \n']);
end;
%format ('1'//10X,'TEST PCH TO B-SPLINE CONVERTER');
if( kprint>=2 )
writef(lun,[ '\n ' , '\n ' ,repmat(' ',1,10),'PCHQK5 RESULTS', '\n ' ,repmat(' ',1,10),'--------------' ' \n']);
end;
%format [/10X,'PCHQK5 RESULTS'/10X,'--------------');
if( kprint>=3 )
writef(lun,[ '\n ' ,repmat(' ',1,4),'(Results should be the same for allmlv KNOTYP values.)' ' \n']);
end;
%format [4X,'(Results should be the same for allmlv KNOTYP values.)');
for knotyp = 2 : -1: -1 ;
[n,x,f,d,dumvar5,knotyp,nknots,t,bcoef,ndim,k,ierr]=pchbs(n,x,f,d,1,knotyp,nknots,t,bcoef,ndim,k,ierr);
if( kprint>=3 )
writef(lun,[ '\n ' ,repmat(' ',1,4),'KNOTYP =','%2i',':  NKNOTS =','%3i',',  NDIM =','%3i',',  K =','%2i',',  IERR =','%3i' ' \n'], knotyp , nknots , ndim , k ,ierr);
end;
%format [4X,'KNOTYP =',i2,':  NKNOTS =',i3,',  NDIM =',i3,',  K =',i2,',  IERR =',i3);
if( ierr~=0 )
ifail = fix(ifail + 1);
if( kprint>=3 )
writef(lun,[' *** Failed -- bad IERR value.' ' \n']);
end;
%format (' *** Failed -- bad IERR value.');
else;
inbv = 1;
fermax = zero;
dermax = zero;
if( kprint>=3 )
writef(lun,[ '\n ' ,repmat(' ',1,15),'X',repmat(' ',1,9),'KNOTS',repmat(' ',1,10),'F',repmat(' ',1,7),'FERR',repmat(' ',1,8),'D',repmat(' ',1,7),'DERR' ' \n']);
%format [15X,'X',9X,'KNOTS',10X,'F',7X,'FERR',8X,'D',7X,'DERR');
writef(lun,[repmat(' ',1,18),repmat('%8.2f',1,2) ' \n'], t(1) , t(2));
j = 1;
end;
for i = 1 : n;
[fcalc ,t,bcoef,ndim,k,dumvar6,x(i),inbv,work]=bvalu(t,bcoef,ndim,k,0,x(i),inbv,work);
ferr = f(i) - fcalc;
relerr= @(err,ansmlv)  abs(err)./max(1.0e-5,abs(ansmlv));
fermax = max(fermax,relerr(ferr,f(i)));
[dcalc ,t,bcoef,ndim,k,dumvar6,x(i),inbv,work]=bvalu(t,bcoef,ndim,k,1,x(i),inbv,work);
derr = d(i) - dcalc;
relerr= @(err,ansmlv)  abs(err)./max(1.0e-5,abs(ansmlv));
dermax = max(dermax,relerr(derr,d(i)));
if( kprint>=3 )
j = fix(j + 2);
writef(lun,[repmat(' ',1,10),repmat('%8.2f',1,3),'%10.4f',repmat('%f',1,1),'%10.2f',repmat('%f',1,0),'%10.4f',repmat('%f',1,1),'%10.2f' ' \n'], x(i) , t(j) , t(j+1) , f(i) , ferr , d(i) ,derr);
%format(10x,3f8.2,f10.4,1p,e10.2,0p,f10.4,1p,e10.2);
end;
end; i = fix(n+1);
if( kprint>=3 )
j = fix(j + 2);
writef(lun,[repmat(' ',1,18),repmat('%8.2f',1,2) ' \n'], t(j) , t(j+1));
end;
fail =(fermax>tol) |(dermax>tol);
if( fail )
ifail = fix(ifail + 1);
end;
if((kprint>=3) ||(kprint>=2) && fail )
writef(lun,[ '\n ' ,repmat(' ',1,5),'Maximum relative errors:', '\n ' ,repmat(' ',1,15),'F-error =',repmat('%f',1,1),'%13.5f',repmat(' ',1,5),'D-error =','%13.5f', '\n ' ,repmat(' ',1,5),'Both should be less than  TOL =','%13.5f' ' \n'],fermax , dermax , tol);
end;
%format [5X,'Maximum relative errors:'/15X,'F-error =',1P,e13.5,5X,'D-error =',e13.5/5X,'Both should be less than  TOL =',e13.5);
end;
if( knotyp==0 )
for i = 1 : nknots;
tsave(i) = t(i);
end; i = fix(nknots+1);
elseif( knotyp==-1 ) ;
termax = zero;
for i = 1 : nknots;
terr = abs(t(i)-tsave(i));
termax = max(termax,terr);
end; i = fix(nknots+1);
if( termax>tolz )
ifail = fix(ifail + 1);
if( kprint>=2 )
writef(lun,[ '\n ' ,' *** T-ARRAY MAXIMUM CHANGE =',repmat('%f',1,1),'%13.5f',';  SHOULD NOT EXCEED TOLZ =','%13.5f' ' \n'], termax , tolz);
end;
%format [' *** T-ARRAY MAXIMUM CHANGE =',1P,e13.5,';  SHOULD NOT EXCEED TOLZ =',e13.5);
end;
end;
end; knotyp = fix(-1 -1);
if((kprint>=2) &&(ifail~=0) )
writef(lun,[ '\n ' ,' *** TROUBLE ***','%5i',' CONVERSION TESTS FAILED.' ' \n'], ifail);
end;
%format [' *** TROUBLE ***',i5,' CONVERSION TESTS FAILED.');
if( ifail==0 )
ipass = 1;
if( kprint>=2 )
writef(lun,[ '\n ' ,' ------------  PCHIP PASSED  ALL CONVERSION TESTS',' ------------' ' \n']);
end;
%format [' ------------  PCHIP PASSED  ALL CONVERSION TESTS',' ------------');
else;
ipass = 0;
if( kprint>=1 )
writef(lun,[ '\n ' ,' **  PCHIP FAILED SOME CONVERSION TESTS',' **' ' \n']);
end;
%format [' **  PCHIP FAILED SOME CONVERSION TESTS',' **');
end;
return;
%format(18x,2f8.2);
end %subroutine pchqk5