function r9atn1 (x)
c april 1978 edition.  w. fullerton c3, los alamos scientific lab.
c
c evaluate  atan(x)  from first order, that is, evaluate
c (atan(x)-x)/x**3  with relative error accuracy so that
c        atan(x) = x + x**3*r9atn1(x).
c
      dimension atn1cs(21)
      external atan, csevl, inits, r1mach, sqrt
c
c series for atn1       on the interval  0.          to  1.00000d+00
c                                        with weighted error   2.21e-17
c                                         log weighted error  16.66
c                               significant figures required  15.44
c                                    decimal places required  17.32
c
      data atn1cs( 1) /   -.0328399753 5355202e0 /
      data atn1cs( 2) /    .0583343234 3172412e0 /
      data atn1cs( 3) /   -.0074003696 9671964e0 /
      data atn1cs( 4) /    .0010097841 9933728e0 /
      data atn1cs( 5) /   -.0001439787 1635652e0 /
      data atn1cs( 6) /    .0000211451 2648992e0 /
      data atn1cs( 7) /   -.0000031723 2107425e0 /
      data atn1cs( 8) /    .0000004836 6203654e0 /
      data atn1cs( 9) /   -.0000000746 7746546e0 /
      data atn1cs(10) /    .0000000116 4800896e0 /
      data atn1cs(11) /   -.0000000018 3208837e0 /
      data atn1cs(12) /    .0000000002 9019082e0 /
      data atn1cs(13) /   -.0000000000 4623885e0 /
      data atn1cs(14) /    .0000000000 0740552e0 /
      data atn1cs(15) /   -.0000000000 0119135e0 /
      data atn1cs(16) /    .0000000000 0019240e0 /
      data atn1cs(17) /   -.0000000000 0003118e0 /
      data atn1cs(18) /    .0000000000 0000506e0 /
      data atn1cs(19) /   -.0000000000 0000082e0 /
      data atn1cs(20) /    .0000000000 0000013e0 /
      data atn1cs(21) /   -.0000000000 0000002e0 /
c
      data ntatn1, xsml, xbig, xmax / 0, 3*0.0 /
c
      if (ntatn1.ne.0) go to 10
      eps = r1mach(3)
      ntatn1 = inits (atn1cs, 21, 0.1*eps)
c
      xsml = sqrt (0.1*eps)
      xbig = 1.571/sqrt(eps)
      xmax = 1.571/eps
c
 10   y = abs(x)
      if (y.gt.1.0) go to 20
c
      if (y.le.xsml) r9atn1 = -1.0/3.0
      if (y.le.xsml) return
c
      r9atn1 = -0.25 + csevl (2.0*y*y-1., atn1cs, ntatn1)
      return
c
 20   if (y.gt.xmax) call seteru (
     1  51hr9atn1  no precision in answer because x is too big, 51,2,2)
      if (y.gt.xbig) call seteru (
     1  53hr9atn1  answer lt half precision because x is too big, 53,
     2  1, 1)
c
      r9atn1 = (atan(x) - x) / x**3
      return
c
      end