function cinh (x)
c december 1980 edition.  w. fullerton, bell labs.
c
c evaluate the hyperbolic cin function.
c         cinh(x) = integral from 0 to x of (cosh(t) - 1)/t dt.
c
      dimension cinhcs(10)
      external alog, chi, csevl, inits, r1mach, sqrt
c
c series for cinh on the interval  0.00000e+00 to  9.00000e+00
c                                        with weighted error   6.64e-17
c                                         log weighted error  16.18
c                               significant figures required  15.08
c                                    decimal places required  16.68
c
      data cinhcs(  1) /     0.1093291636 520734431 e0/
      data cinhcs(  2) /     0.0573928847 550379676 e0/
      data cinhcs(  3) /     0.0028095756 978830353 e0/
      data cinhcs(  4) /     0.0000828780 840721357 e0/
      data cinhcs(  5) /     0.0000016278 596173914 e0/
      data cinhcs(  6) /     0.0000000227 809519256 e0/
      data cinhcs(  7) /     0.0000000002 384484842 e0/
      data cinhcs(  8) /     0.0000000000 019360830 e0/
      data cinhcs(  9) /     0.0000000000 000125454 e0/
      data cinhcs( 10) /     0.0000000000 000000664 e0/
c
      data eul / 0.5772156649 0153286 e0 /
      data ncinh, xsml, xmin /0, 2*0.0/
c
      if (ncinh.ne.0) go to 10
      ncinh = inits (cinhcs, 10, 0.1*r1mach(3))
      xsml = sqrt (r1mach(3))
      xmin = 2.0*sqrt(r1mach(1))
c
 10   absx = abs(x)
      if (absx.gt.3.0) go to 20
c
      cinh = 0.0
      if (x.ne.0.0 .and. absx.le.xmin) call seteru (
     1  39hcinh    abs(x) so small cinh underflows, 39, 1, 0)
      if (absx.le.xmin) return
c
      y = -1.0
      if (absx.gt.xsml) y = x*x/9.0 - 1.0
      cinh = x*x* (0.25 + csevl (y, cinhcs, ncinh))
      return
c
 20   cinh = chi(absx) - eul - alog(absx)
c
      return
      end