% Fully developed turbulent flow in smooth wall pipe
function [h, NuD, f, ReD] = TurbPipeNu(U,D,mat,Tm)
  [ReD Pr km] = Re(U,D,mat,Tm);                                 
  if     Pr < 0.5    warning = 'TurbPipeNu: Pr < 0.5'
  elseif Pr > 2000   warning = 'TurbPipeNu: Pr > 2000' 
  end
  % published limit of ReD > 3000 has been changed to 2300 below
  if min(ReD) < 2300 warning = 'TurbPipeNu: ReD < 2300' 
  end
  f   = (0.8686*log( ReD./(1.964*log(ReD)-3.8215) )).^-2;       
  NuD = (f/8).*(ReD-1000).*Pr./(1+12.7*sqrt(f/8)*(Pr^(2/3)-1)); 
  NuD(find(NuD<1))=1;                                           % (in case ReD<1000)
  h   = NuD*km./D;
end
