src/numerical/slatec/zseri.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) (coner 1.0) (conei 0.0))
  21   (declare (type (double-float) zeror zeroi coner conei))
  22   (defun zseri (zr zi fnu kode n yr yi nz tol elim alim)
  23     (declare (type (simple-array double-float (*)) yi yr)
  24              (type (f2cl-lib:integer4) nz n kode)
  25              (type (double-float) alim elim tol fnu zi zr))
  26     (prog ((wr (make-array 2 :element-type 'double-float))
  27            (wi (make-array 2 :element-type 'double-float)) (i 0) (ib 0)
  28            (idum 0) (iflag 0) (il 0) (k 0) (l 0) (m 0) (nn 0) (nw 0) (aa 0.0)
  29            (acz 0.0) (ak 0.0) (ak1i 0.0) (ak1r 0.0) (arm 0.0) (ascle 0.0)
  30            (atol 0.0) (az 0.0) (cki 0.0) (ckr 0.0) (coefi 0.0) (coefr 0.0)
  31            (crscr 0.0) (czi 0.0) (czr 0.0) (dfnu 0.0) (fnup 0.0) (hzi 0.0)
  32            (hzr 0.0) (raz 0.0) (rs 0.0) (rtr1 0.0) (rzi 0.0) (rzr 0.0) (s 0.0)
  33            (ss 0.0) (sti 0.0) (str 0.0) (s1i 0.0) (s1r 0.0) (s2i 0.0)
  34            (s2r 0.0))
  35       (declare (type (simple-array double-float (2)) wr wi)
  36                (type (double-float) s2r s2i s1r s1i str sti ss s rzr rzi rtr1
  37                                     rs raz hzr hzi fnup dfnu czr czi crscr
  38                                     coefr coefi ckr cki az atol ascle arm ak1r
  39                                     ak1i ak acz aa)
  40                (type (f2cl-lib:integer4) nw nn m l k il iflag idum ib i))
  41       (setf nz 0)
  42       (setf az (coerce (realpart (zabs zr zi)) 'double-float))
  43       (if (= az 0.0) (go label160))
  44       (setf arm (* 1000.0 (f2cl-lib:d1mach 1)))
  45       (setf rtr1 (f2cl-lib:fsqrt arm))
  46       (setf crscr 1.0)
  47       (setf iflag 0)
  48       (if (< az arm) (go label150))
  49       (setf hzr (* 0.5 zr))
  50       (setf hzi (* 0.5 zi))
  51       (setf czr zeror)
  52       (setf czi zeroi)
  53       (if (<= az rtr1) (go label10))
  54       (multiple-value-bind (var-0 var-1 var-2 var-3 var-4 var-5)
  55           (zmlt hzr hzi hzr hzi czr czi)
  56         (declare (ignore var-0 var-1 var-2 var-3))
  57         (setf czr var-4)
  58         (setf czi var-5))
  59      label10
  60       (setf acz (coerce (realpart (zabs czr czi)) 'double-float))
  61       (setf nn n)
  62       (multiple-value-bind (var-0 var-1 var-2 var-3 var-4)
  63           (zlog hzr hzi ckr cki idum)
  64         (declare (ignore var-0 var-1))
  65         (setf ckr var-2)
  66         (setf cki var-3)
  67         (setf idum var-4))
  68      label20
  69       (setf dfnu (+ fnu (f2cl-lib:int-sub nn 1)))
  70       (setf fnup (+ dfnu 1.0))
  71       (setf ak1r (* ckr dfnu))
  72       (setf ak1i (* cki dfnu))
  73       (setf ak
  74               (multiple-value-bind (ret-val var-0 var-1)
  75                   (dgamln fnup idum)
  76                 (declare (ignore var-0))
  77                 (setf idum var-1)
  78                 ret-val))
  79       (setf ak1r (- ak1r ak))
  80       (if (= kode 2) (setf ak1r (- ak1r zr)))
  81       (if (> ak1r (- elim)) (go label40))
  82      label30
  83       (setf nz (f2cl-lib:int-add nz 1))
  84       (setf (f2cl-lib:fref yr (nn) ((1 n))) zeror)
  85       (setf (f2cl-lib:fref yi (nn) ((1 n))) zeroi)
  86       (if (> acz dfnu) (go label190))
  87       (setf nn (f2cl-lib:int-sub nn 1))
  88       (if (= nn 0) (go end_label))
  89       (go label20)
  90      label40
  91       (if (> ak1r (- alim)) (go label50))
  92       (setf iflag 1)
  93       (setf ss (/ 1.0 tol))
  94       (setf crscr tol)
  95       (setf ascle (* arm ss))
  96      label50
  97       (setf aa (exp ak1r))
  98       (if (= iflag 1) (setf aa (* aa ss)))
  99       (setf coefr (* aa (cos ak1i)))
 100       (setf coefi (* aa (sin ak1i)))
 101       (setf atol (/ (* tol acz) fnup))
 102       (setf il (min (the f2cl-lib:integer4 2) (the f2cl-lib:integer4 nn)))
 103       (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1))
 104                     ((> i il) nil)
 105         (tagbody
 106           (setf dfnu (+ fnu (f2cl-lib:int-sub nn i)))
 107           (setf fnup (+ dfnu 1.0))
 108           (setf s1r coner)
 109           (setf s1i conei)
 110           (if (< acz (* tol fnup)) (go label70))
 111           (setf ak1r coner)
 112           (setf ak1i conei)
 113           (setf ak (+ fnup 2.0))
 114           (setf s fnup)
 115           (setf aa 2.0)
 116          label60
 117           (setf rs (/ 1.0 s))
 118           (setf str (- (* ak1r czr) (* ak1i czi)))
 119           (setf sti (+ (* ak1r czi) (* ak1i czr)))
 120           (setf ak1r (* str rs))
 121           (setf ak1i (* sti rs))
 122           (setf s1r (+ s1r ak1r))
 123           (setf s1i (+ s1i ak1i))
 124           (setf s (+ s ak))
 125           (setf ak (+ ak 2.0))
 126           (setf aa (* aa acz rs))
 127           (if (> aa atol) (go label60))
 128          label70
 129           (setf s2r (- (* s1r coefr) (* s1i coefi)))
 130           (setf s2i (+ (* s1r coefi) (* s1i coefr)))
 131           (setf (f2cl-lib:fref wr (i) ((1 2))) s2r)
 132           (setf (f2cl-lib:fref wi (i) ((1 2))) s2i)
 133           (if (= iflag 0) (go label80))
 134           (multiple-value-bind (var-0 var-1 var-2 var-3 var-4)
 135               (zuchk s2r s2i nw ascle tol)
 136             (declare (ignore var-0 var-1 var-3 var-4))
 137             (setf nw var-2))
 138           (if (/= nw 0) (go label30))
 139          label80
 140           (setf m (f2cl-lib:int-add (f2cl-lib:int-sub nn i) 1))
 141           (setf (f2cl-lib:fref yr (m) ((1 n))) (* s2r crscr))
 142           (setf (f2cl-lib:fref yi (m) ((1 n))) (* s2i crscr))
 143           (if (= i il) (go label90))
 144           (multiple-value-bind (var-0 var-1 var-2 var-3 var-4 var-5)
 145               (zdiv coefr coefi hzr hzi str sti)
 146             (declare (ignore var-0 var-1 var-2 var-3))
 147             (setf str var-4)
 148             (setf sti var-5))
 149           (setf coefr (* str dfnu))
 150           (setf coefi (* sti dfnu))
 151          label90))
 152       (if (<= nn 2) (go end_label))
 153       (setf k (f2cl-lib:int-sub nn 2))
 154       (setf ak (coerce (the f2cl-lib:integer4 k) 'double-float))
 155       (setf raz (/ 1.0 az))
 156       (setf str (* zr raz))
 157       (setf sti (* (- zi) raz))
 158       (setf rzr (* (+ str str) raz))
 159       (setf rzi (* (+ sti sti) raz))
 160       (if (= iflag 1) (go label120))
 161       (setf ib 3)
 162      label100
 163       (f2cl-lib:fdo (i ib (f2cl-lib:int-add i 1))
 164                     ((> i nn) nil)
 165         (tagbody
 166           (setf (f2cl-lib:fref yr (k) ((1 n)))
 167                   (+
 168                    (* (+ ak fnu)
 169                       (-
 170                        (* rzr
 171                           (f2cl-lib:fref yr ((f2cl-lib:int-add k 1)) ((1 n))))
 172                        (* rzi
 173                           (f2cl-lib:fref yi
 174                                          ((f2cl-lib:int-add k 1))
 175                                          ((1 n))))))
 176                    (f2cl-lib:fref yr ((f2cl-lib:int-add k 2)) ((1 n)))))
 177           (setf (f2cl-lib:fref yi (k) ((1 n)))
 178                   (+
 179                    (* (+ ak fnu)
 180                       (+
 181                        (* rzr
 182                           (f2cl-lib:fref yi ((f2cl-lib:int-add k 1)) ((1 n))))
 183                        (* rzi
 184                           (f2cl-lib:fref yr
 185                                          ((f2cl-lib:int-add k 1))
 186                                          ((1 n))))))
 187                    (f2cl-lib:fref yi ((f2cl-lib:int-add k 2)) ((1 n)))))
 188           (setf ak (- ak 1.0))
 189           (setf k (f2cl-lib:int-sub k 1))
 190          label110))
 191       (go end_label)
 192      label120
 193       (setf s1r (f2cl-lib:fref wr (1) ((1 2))))
 194       (setf s1i (f2cl-lib:fref wi (1) ((1 2))))
 195       (setf s2r (f2cl-lib:fref wr (2) ((1 2))))
 196       (setf s2i (f2cl-lib:fref wi (2) ((1 2))))
 197       (f2cl-lib:fdo (l 3 (f2cl-lib:int-add l 1))
 198                     ((> l nn) nil)
 199         (tagbody
 200           (setf ckr s2r)
 201           (setf cki s2i)
 202           (setf s2r (+ s1r (* (+ ak fnu) (- (* rzr ckr) (* rzi cki)))))
 203           (setf s2i (+ s1i (* (+ ak fnu) (+ (* rzr cki) (* rzi ckr)))))
 204           (setf s1r ckr)
 205           (setf s1i cki)
 206           (setf ckr (* s2r crscr))
 207           (setf cki (* s2i crscr))
 208           (setf (f2cl-lib:fref yr (k) ((1 n))) ckr)
 209           (setf (f2cl-lib:fref yi (k) ((1 n))) cki)
 210           (setf ak (- ak 1.0))
 211           (setf k (f2cl-lib:int-sub k 1))
 212           (if (> (zabs ckr cki) ascle) (go label140))
 213          label130))
 214       (go end_label)
 215      label140
 216       (setf ib (f2cl-lib:int-add l 1))
 217       (if (> ib nn) (go end_label))
 218       (go label100)
 219      label150
 220       (setf nz n)
 221       (if (= fnu 0.0) (setf nz (f2cl-lib:int-sub nz 1)))
 222      label160
 223       (setf (f2cl-lib:fref yr (1) ((1 n))) zeror)
 224       (setf (f2cl-lib:fref yi (1) ((1 n))) zeroi)
 225       (if (/= fnu 0.0) (go label170))
 226       (setf (f2cl-lib:fref yr (1) ((1 n))) coner)
 227       (setf (f2cl-lib:fref yi (1) ((1 n))) conei)
 228      label170
 229       (if (= n 1) (go end_label))
 230       (f2cl-lib:fdo (i 2 (f2cl-lib:int-add i 1))
 231                     ((> i n) nil)
 232         (tagbody
 233           (setf (f2cl-lib:fref yr (i) ((1 n))) zeror)
 234           (setf (f2cl-lib:fref yi (i) ((1 n))) zeroi)
 235          label180))
 236       (go end_label)
 237      label190
 238       (setf nz (f2cl-lib:int-sub nz))
 239       (go end_label)
 240      end_label
 241       (return (values nil nil nil nil nil nil nil nz nil nil nil)))))
 242
 243 (in-package #:cl-user)
 244 #+#.(cl:if (cl:find-package '#:f2cl) '(and) '(or))
 245 (eval-when (:load-toplevel :compile-toplevel :execute)
 246   (setf (gethash 'fortran-to-lisp::zseri fortran-to-lisp::*f2cl-function-info*)
 247           (fortran-to-lisp::make-f2cl-finfo
 248            :arg-types '((double-float) (double-float) (double-float)
 249                         (fortran-to-lisp::integer4) (fortran-to-lisp::integer4)
 250                         (simple-array double-float (*))
 251                         (simple-array double-float (*))
 252                         (fortran-to-lisp::integer4) (double-float)
 253                         (double-float) (double-float))
 254            :return-values '(nil nil nil nil nil nil nil fortran-to-lisp::nz nil
 255                             nil nil)
 256            :calls '(fortran-to-lisp::zdiv fortran-to-lisp::zuchk
 257                     fortran-to-lisp::dgamln fortran-to-lisp::zlog
 258                     fortran-to-lisp::zmlt fortran-to-lisp::d1mach
 259                     fortran-to-lisp::zabs))))
 260