function [ab,dig]=dig_Enufinmod(N,nu,c,dig0,dd,nofdig)
%DIG_ENUFINMOD Number of digits needed to obtain polynomials
%orthogonal with respect to the subrange exponential integral.
%   [AB,DIG]=DIG_ENUFINMOD(N,NU,C,DIG0,DD,NOFDIG) determines the 
%   number DIG of digits needed to obtain the first N recurrence
%   coefficients for the weight function w(x)=E_NU(x) on (0,C], C>0, 
%   to an accuracy of NOFDIG decimal digits, using the modified
%   Chebyshev algorithm relative to monic shifted Jacobi polynomials 
%   with Jacobi parameters alpha=0, beta=NU-1 transformed from the
%   interval [0,1] to the interval [0,C]. The routine successively 
%   increments an initial estimate DIG0 by DD units until the
%   required accuracy is achieved. The recurrence coefficients to
%   NOFDIG digits are placed into the Nx2 array AB.

tic
disp('Be sure to select the correct mom0.')
mom0=loadvpa('mmom_Enufin_1_1.txt',9,50,80,1);
%mom0=loadvpa('mmom_Enufin_1_3.txt',9,190,80,1); 
%mom0=loadvpa('mmom_Enufin_1_5.txt',9,190,80,1); 
%mom0=loadvpa('mmom_Enufin_1_10.txt',9,190,80,1);
%mom0=loadvpa('mmom_Enufin_1_20.txt',9,190,80,1);
i=dig0-dd; 
maxerr_Enufin=1;
while maxerr_Enufin>.5*10^(-nofdig)
  i=i+dd; dig=i
  digits(dig); nu=vpa(nu); c=vpa(c);
  abmm0=sr_jacobi01(dig,2*N-1,0,nu-1);
  abmm(:,1)=c*abmm0(:,1); abmm(1,2)=(c^nu)*abmm0(1,2);
  abmm(2:2*N-1,2)=c^2*abmm0(2:2*N-1,2); 
  mom=vpa(mom0,dig);
  if i==dig0
    ab0=schebyshev(dig,N,mom',abmm);
  else
    ab1=schebyshev(dig,N,mom',abmm);
    err=abs((ab1-ab0)./ab1); err=double(err);
    maxerr_Enufin=max(max(err)) 
    ab0=ab1;
  end 
end
ab=vpa(ab1,nofdig);
toc
