function ab=sr_milovanovic(dig,N,a) 
%Recurrence coefficients for the Milovanovic weight function.
%   AB=SR_MILOVANOVIC(DIG,N,A) generates in DIG-digit
%   arithmetic the Nx2 array AB of the first N recurrence
%   coefficients for the weight function x(x)=exp(-x^(-a)-x^a)
%   on [0,Inf], A>0.

digits(dig); a=vpa(a);
mom=vpa(zeros(1,2*N)); 
for k=1:2*N
  mom(k)=vpa(2/a)*besselk(vpa(k/a),2);
end
mom'
ab=schebyshev(dig,N,mom)
