docs: def _gser(a,x):
      "Series representation of Gamma. NumRec sect 6.1."
      ITMAX=100
      EPS=3.e-7
  
      gln=lgamma(a)
      assert(x>=0),'x < 0 in gser'
      if x == 0 : return 0,gln
  
      ap = a
      delt = sum = 1./a
      for i in range(ITMAX):
          ap=ap+1.
          delt=delt*x/ap
          sum=sum+delt
          if abs(delt) < abs(sum)*EPS: break
      else:
          print('a too large, ITMAX too small in gser')
      gamser=sum*np.exp(-x+a*np.log(x)-gln)
      return gamser,gln