src/numerical/slatec/zkscl.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 ((zeror 0.0) (zeroi 0.0))
  21   (declare (type (double-float) zeror zeroi))
  22   (defun zkscl (zrr zri fnu n yr yi nz rzr rzi ascle tol elim)
  23     (declare (type (simple-array double-float (*)) yi yr)
  24              (type (f2cl-lib:integer4) nz n)
  25              (type (double-float) elim tol ascle rzi rzr fnu zri zrr))
  26     (prog ((cyr (make-array 2 :element-type 'double-float))
  27            (cyi (make-array 2 :element-type 'double-float)) (i 0) (ic 0)
  28            (idum 0) (kk 0) (nn 0) (nw 0) (acs 0.0) (as 0.0) (cki 0.0) (ckr 0.0)
  29            (csi 0.0) (csr 0.0) (fn 0.0) (str 0.0) (s1i 0.0) (s1r 0.0) (s2i 0.0)
  30            (s2r 0.0) (zdr 0.0) (zdi 0.0) (celmr 0.0) (elm 0.0) (helim 0.0)
  31            (alas 0.0))
  32       (declare (type (simple-array double-float (2)) cyr cyi)
  33                (type (double-float) alas helim elm celmr zdi zdr s2r s2i s1r
  34                                     s1i str fn csr csi ckr cki as acs)
  35                (type (f2cl-lib:integer4) nw nn kk idum ic i))
  36       (setf nz 0)
  37       (setf ic 0)
  38       (setf nn (min (the f2cl-lib:integer4 2) (the f2cl-lib:integer4 n)))
  39       (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1))
  40                     ((> i nn) nil)
  41         (tagbody
  42           (setf s1r (f2cl-lib:fref yr (i) ((1 n))))
  43           (setf s1i (f2cl-lib:fref yi (i) ((1 n))))
  44           (setf (f2cl-lib:fref cyr (i) ((1 2))) s1r)
  45           (setf (f2cl-lib:fref cyi (i) ((1 2))) s1i)
  46           (setf as (coerce (realpart (zabs s1r s1i)) 'double-float))
  47           (setf acs (- (f2cl-lib:flog as) zrr))
  48           (setf nz (f2cl-lib:int-add nz 1))
  49           (setf (f2cl-lib:fref yr (i) ((1 n))) zeror)
  50           (setf (f2cl-lib:fref yi (i) ((1 n))) zeroi)
  51           (if (< acs (- elim)) (go label10))
  52           (multiple-value-bind (var-0 var-1 var-2 var-3 var-4)
  53               (zlog s1r s1i csr csi idum)
  54             (declare (ignore var-0 var-1))
  55             (setf csr var-2)
  56             (setf csi var-3)
  57             (setf idum var-4))
  58           (setf csr (- csr zrr))
  59           (setf csi (- csi zri))
  60           (setf str (/ (exp csr) tol))
  61           (setf csr (* str (cos csi)))
  62           (setf csi (* str (sin csi)))
  63           (multiple-value-bind (var-0 var-1 var-2 var-3 var-4)
  64               (zuchk csr csi nw ascle tol)
  65             (declare (ignore var-0 var-1 var-3 var-4))
  66             (setf nw var-2))
  67           (if (/= nw 0) (go label10))
  68           (setf (f2cl-lib:fref yr (i) ((1 n))) csr)
  69           (setf (f2cl-lib:fref yi (i) ((1 n))) csi)
  70           (setf ic i)
  71           (setf nz (f2cl-lib:int-sub nz 1))
  72          label10))
  73       (if (= n 1) (go end_label))
  74       (if (> ic 1) (go label20))
  75       (setf (f2cl-lib:fref yr (1) ((1 n))) zeror)
  76       (setf (f2cl-lib:fref yi (1) ((1 n))) zeroi)
  77       (setf nz 2)
  78      label20
  79       (if (= n 2) (go end_label))
  80       (if (= nz 0) (go end_label))
  81       (setf fn (+ fnu 1.0))
  82       (setf ckr (* fn rzr))
  83       (setf cki (* fn rzi))
  84       (setf s1r (f2cl-lib:fref cyr (1) ((1 2))))
  85       (setf s1i (f2cl-lib:fref cyi (1) ((1 2))))
  86       (setf s2r (f2cl-lib:fref cyr (2) ((1 2))))
  87       (setf s2i (f2cl-lib:fref cyi (2) ((1 2))))
  88       (setf helim (* 0.5 elim))
  89       (setf elm (exp (- elim)))
  90       (setf celmr elm)
  91       (setf zdr zrr)
  92       (setf zdi zri)
  93       (f2cl-lib:fdo (i 3 (f2cl-lib:int-add i 1))
  94                     ((> i n) nil)
  95         (tagbody
  96           (setf kk i)
  97           (setf csr s2r)
  98           (setf csi s2i)
  99           (setf s2r (+ (- (* ckr csr) (* cki csi)) s1r))
 100           (setf s2i (+ (* cki csr) (* ckr csi) s1i))
 101           (setf s1r csr)
 102           (setf s1i csi)
 103           (setf ckr (+ ckr rzr))
 104           (setf cki (+ cki rzi))
 105           (setf as (coerce (realpart (zabs s2r s2i)) 'double-float))
 106           (setf alas (f2cl-lib:flog as))
 107           (setf acs (- alas zdr))
 108           (setf nz (f2cl-lib:int-add nz 1))
 109           (setf (f2cl-lib:fref yr (i) ((1 n))) zeror)
 110           (setf (f2cl-lib:fref yi (i) ((1 n))) zeroi)
 111           (if (< acs (- elim)) (go label25))
 112           (multiple-value-bind (var-0 var-1 var-2 var-3 var-4)
 113               (zlog s2r s2i csr csi idum)
 114             (declare (ignore var-0 var-1))
 115             (setf csr var-2)
 116             (setf csi var-3)
 117             (setf idum var-4))
 118           (setf csr (- csr zdr))
 119           (setf csi (- csi zdi))
 120           (setf str (/ (exp csr) tol))
 121           (setf csr (* str (cos csi)))
 122           (setf csi (* str (sin csi)))
 123           (multiple-value-bind (var-0 var-1 var-2 var-3 var-4)
 124               (zuchk csr csi nw ascle tol)
 125             (declare (ignore var-0 var-1 var-3 var-4))
 126             (setf nw var-2))
 127           (if (/= nw 0) (go label25))
 128           (setf (f2cl-lib:fref yr (i) ((1 n))) csr)
 129           (setf (f2cl-lib:fref yi (i) ((1 n))) csi)
 130           (setf nz (f2cl-lib:int-sub nz 1))
 131           (if (= ic (f2cl-lib:int-sub kk 1)) (go label40))
 132           (setf ic kk)
 133           (go label30)
 134          label25
 135           (if (< alas helim) (go label30))
 136           (setf zdr (- zdr elim))
 137           (setf s1r (* s1r celmr))
 138           (setf s1i (* s1i celmr))
 139           (setf s2r (* s2r celmr))
 140           (setf s2i (* s2i celmr))
 141          label30))
 142       (setf nz n)
 143       (if (= ic n) (setf nz (f2cl-lib:int-sub n 1)))
 144       (go label45)
 145      label40
 146       (setf nz (f2cl-lib:int-sub kk 2))
 147      label45
 148       (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1))
 149                     ((> i nz) nil)
 150         (tagbody
 151           (setf (f2cl-lib:fref yr (i) ((1 n))) zeror)
 152           (setf (f2cl-lib:fref yi (i) ((1 n))) zeroi)
 153          label50))
 154       (go end_label)
 155      end_label
 156       (return (values nil nil nil nil nil nil nz nil nil nil nil nil)))))
 157
 158 (in-package #:cl-user)
 159 #+#.(cl:if (cl:find-package '#:f2cl) '(and) '(or))
 160 (eval-when (:load-toplevel :compile-toplevel :execute)
 161   (setf (gethash 'fortran-to-lisp::zkscl fortran-to-lisp::*f2cl-function-info*)
 162           (fortran-to-lisp::make-f2cl-finfo
 163            :arg-types '((double-float) (double-float) (double-float)
 164                         (fortran-to-lisp::integer4)
 165                         (simple-array double-float (*))
 166                         (simple-array double-float (*))
 167                         (fortran-to-lisp::integer4) (double-float)
 168                         (double-float) (double-float) (double-float)
 169                         (double-float))
 170            :return-values '(nil nil nil nil nil nil fortran-to-lisp::nz nil nil
 171                             nil nil nil)
 172            :calls '(fortran-to-lisp::zuchk fortran-to-lisp::zlog
 173                     fortran-to-lisp::zabs))))
 174