src/numerical/slatec/dlngam.lisp

   1 ;;; Compiled by f2cl version:
   2 ;;; ("f2cl1.l,v 46c1f6a93b0d 2012/05/03 04:40:28 toy $"
   3 ;;;  "f2cl2.l,v 96616d88fb7e 2008/02/22 22:19:34 rtoy $"
   4 ;;;  "f2cl3.l,v 96616d88fb7e 2008/02/22 22:19:34 rtoy $"
   5 ;;;  "f2cl4.l,v 96616d88fb7e 2008/02/22 22:19:34 rtoy $"
   6 ;;;  "f2cl5.l,v 46c1f6a93b0d 2012/05/03 04:40:28 toy $"
   7 ;;;  "f2cl6.l,v 1d5cbacbb977 2008/08/24 00:56:27 rtoy $"
   8 ;;;  "macros.l,v fceac530ef0c 2011/11/26 04:02:26 toy $")
   9
  10 ;;; Using Lisp CMU Common Lisp snapshot-2012-04 (20C Unicode)
  11 ;;;
  12 ;;; Options: ((:prune-labels nil) (:auto-save t) (:relaxed-array-decls t)
  13 ;;;           (:coerce-assigns :as-needed) (:array-type ':simple-array)
  14 ;;;           (:array-slicing nil) (:declare-common nil)
  15 ;;;           (:float-format double-float))
  16
  17 (in-package :slatec)
  18
  19
  20 (let ((xmax 0.0)
  21       (dxrel 0.0)
  22       (sq2pil 0.9189385332046728)
  23       (sqpi2l 0.22579135264472744)
  24       (pi$ 3.141592653589793)
  25       (first$ nil))
  26   (declare (type (double-float) xmax dxrel sq2pil sqpi2l pi$)
  27            (type f2cl-lib:logical first$))
  28   (setq first$ f2cl-lib:%true%)
  29   (defun dlngam (x)
  30     (declare (type (double-float) x))
  31     (prog ((sinpiy 0.0) (y 0.0) (temp 0.0) (dlngam 0.0))
  32       (declare (type (double-float) dlngam temp y sinpiy))
  33       (cond
  34         (first$
  35          (setf temp (/ 1.0 (f2cl-lib:flog (f2cl-lib:d1mach 2))))
  36          (setf xmax (* temp (f2cl-lib:d1mach 2)))
  37          (setf dxrel (f2cl-lib:fsqrt (f2cl-lib:d1mach 4)))))
  38       (setf first$ f2cl-lib:%false%)
  39       (setf y (abs x))
  40       (if (> y 10.0) (go label20))
  41       (setf dlngam (f2cl-lib:flog (abs (dgamma x))))
  42       (go end_label)
  43      label20
  44       (if (> y xmax)
  45           (xermsg "SLATEC" "DLNGAM" "ABS(X) SO BIG DLNGAM OVERFLOWS" 2 2))
  46       (if (> x 0.0)
  47           (setf dlngam
  48                   (+ (- (+ sq2pil (* (- x 0.5) (f2cl-lib:flog x))) x)
  49                      (d9lgmc y))))
  50       (if (> x 0.0) (go end_label))
  51       (setf sinpiy (abs (sin (* pi$ y))))
  52       (if (= sinpiy 0.0)
  53           (xermsg "SLATEC" "DLNGAM" "X IS A NEGATIVE INTEGER" 3 2))
  54       (if (< (abs (/ (- x (f2cl-lib:aint (- x 0.5))) x)) dxrel)
  55           (xermsg "SLATEC" "DLNGAM"
  56            "ANSWER LT HALF PRECISION BECAUSE X TOO NEAR NEGATIVE INTEGER" 1 1))
  57       (setf dlngam
  58               (- (+ sqpi2l (* (- x 0.5) (f2cl-lib:flog y)))
  59                  x
  60                  (f2cl-lib:flog sinpiy)
  61                  (d9lgmc y)))
  62       (go end_label)
  63      end_label
  64       (return (values dlngam nil)))))
  65
  66 (in-package #:cl-user)
  67 #+#.(cl:if (cl:find-package '#:f2cl) '(and) '(or))
  68 (eval-when (:load-toplevel :compile-toplevel :execute)
  69   (setf (gethash 'fortran-to-lisp::dlngam
  70                  fortran-to-lisp::*f2cl-function-info*)
  71           (fortran-to-lisp::make-f2cl-finfo :arg-types '((double-float))
  72                                             :return-values '(nil)
  73                                             :calls '(fortran-to-lisp::d9lgmc
  74                                                      fortran-to-lisp::xermsg
  75                                                      fortran-to-lisp::dgamma
  76                                                      fortran-to-lisp::d1mach))))
  77