share/odepack/src/dpcg.lisp

   1 ;;; Compiled by f2cl version:
   2 ;;; ("f2cl1.l,v 95098eb54f13 2013/04/01 00:45:16 toy $"
   3 ;;;  "f2cl2.l,v 95098eb54f13 2013/04/01 00:45:16 toy $"
   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 95098eb54f13 2013/04/01 00:45:16 toy $"
   7 ;;;  "f2cl6.l,v 1d5cbacbb977 2008/08/24 00:56:27 rtoy $"
   8 ;;;  "macros.l,v 1409c1352feb 2013/03/24 20:44:50 toy $")
   9
  10 ;;; Using Lisp CMU Common Lisp snapshot-2013-11 (20E Unicode)
  11 ;;;
  12 ;;; Options: ((:prune-labels nil) (:auto-save t) (:relaxed-array-decls t)
  13 ;;;           (:coerce-assigns :as-needed) (:array-type ':array)
  14 ;;;           (:array-slicing t) (:declare-common nil)
  15 ;;;           (:float-format single-float))
  16
  17 (in-package "ODEPACK")
  18
  19
  20 (defun dpcg
  21        (neq tn y savf r wght n maxl delta hl0 jpre mnewt f psol npsl x p w z
  22         lpcg wp iwp wk iflag)
  23   (declare (type (f2cl-lib:integer4) iflag lpcg npsl mnewt jpre maxl n)
  24            (type (array double-float (*)) wk wp z w p x wght r savf y)
  25            (type (double-float) hl0 delta tn)
  26            (type (array f2cl-lib:integer4 (*)) iwp neq))
  27   (f2cl-lib:with-multi-array-data
  28       ((neq f2cl-lib:integer4 neq-%data% neq-%offset%)
  29        (iwp f2cl-lib:integer4 iwp-%data% iwp-%offset%)
  30        (y double-float y-%data% y-%offset%)
  31        (savf double-float savf-%data% savf-%offset%)
  32        (r double-float r-%data% r-%offset%)
  33        (wght double-float wght-%data% wght-%offset%)
  34        (x double-float x-%data% x-%offset%)
  35        (p double-float p-%data% p-%offset%)
  36        (w double-float w-%data% w-%offset%)
  37        (z double-float z-%data% z-%offset%)
  38        (wp double-float wp-%data% wp-%offset%)
  39        (wk double-float wk-%data% wk-%offset%))
  40     (prog ((alpha 0.0d0) (beta 0.0d0) (bnrm 0.0d0) (ptw 0.0d0) (rnrm 0.0d0)
  41            (ztr 0.0d0) (ztr0 0.0d0) (i 0) (ier 0))
  42       (declare (type (f2cl-lib:integer4) ier i)
  43                (type (double-float) ztr0 ztr rnrm ptw bnrm beta alpha))
  44       (setf iflag 0)
  45       (setf npsl 0)
  46       (setf lpcg 0)
  47       (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1))
  48                     ((> i n) nil)
  49         (tagbody
  50          label10
  51           (setf (f2cl-lib:fref x-%data% (i) ((1 *)) x-%offset%) 0.0d0)))
  52       (setf bnrm (dvnorm n r wght))
  53       (if (> bnrm delta) (go label20))
  54       (if (> mnewt 0) (go end_label))
  55       (dcopy n r 1 x 1)
  56       (go end_label)
  57      label20
  58       (setf ztr 0.0d0)
  59      label30
  60       (setf lpcg (f2cl-lib:int-add lpcg 1))
  61       (dcopy n r 1 z 1)
  62       (setf ier 0)
  63       (if (= jpre 0) (go label40))
  64       (multiple-value-bind
  65             (var-0 var-1 var-2 var-3 var-4 var-5 var-6 var-7 var-8 var-9
  66              var-10)
  67           (funcall psol neq tn y savf wk hl0 wp iwp z 3 ier)
  68         (declare (ignore var-0 var-2 var-3 var-4 var-6 var-7 var-8 var-9))
  69         (when var-1
  70           (setf tn var-1))
  71         (when var-5
  72           (setf hl0 var-5))
  73         (when var-10
  74           (setf ier var-10)))
  75       (setf npsl (f2cl-lib:int-add npsl 1))
  76       (if (/= ier 0) (go label100))
  77      label40
  78       (setf ztr0 ztr)
  79       (setf ztr (ddot n z 1 r 1))
  80       (if (/= lpcg 1) (go label50))
  81       (dcopy n z 1 p 1)
  82       (go label70)
  83      label50
  84       (if (= ztr0 0.0d0) (go label200))
  85       (setf beta (/ ztr ztr0))
  86       (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1))
  87                     ((> i n) nil)
  88         (tagbody
  89          label60
  90           (setf (f2cl-lib:fref p-%data% (i) ((1 *)) p-%offset%)
  91                   (+ (f2cl-lib:fref z-%data% (i) ((1 *)) z-%offset%)
  92                      (* beta
  93                         (f2cl-lib:fref p-%data% (i) ((1 *)) p-%offset%))))))
  94      label70
  95       (datp neq y savf p wght hl0 wk f w)
  96       (setf ptw (ddot n p 1 w 1))
  97       (if (= ptw 0.0d0) (go label200))
  98       (setf alpha (/ ztr ptw))
  99       (daxpy n alpha p 1 x 1)
 100       (setf alpha (- alpha))
 101       (daxpy n alpha w 1 r 1)
 102       (setf rnrm (dvnorm n r wght))
 103       (if (<= rnrm delta) (go end_label))
 104       (if (< lpcg maxl) (go label30))
 105       (setf iflag 2)
 106       (if (<= rnrm 1.0d0) (setf iflag 1))
 107       (if (and (<= rnrm bnrm) (= mnewt 0)) (setf iflag 1))
 108       (go end_label)
 109      label100
 110       (if (< ier 0) (setf iflag -1))
 111       (if (> ier 0) (setf iflag 3))
 112       (go end_label)
 113      label200
 114       (setf iflag 4)
 115       (go end_label)
 116      end_label
 117       (return
 118        (values nil
 119                tn
 120                nil
 121                nil
 122                nil
 123                nil
 124                nil
 125                nil
 126                nil
 127                hl0
 128                nil
 129                nil
 130                nil
 131                nil
 132                npsl
 133                nil
 134                nil
 135                nil
 136                nil
 137                lpcg
 138                nil
 139                nil
 140                nil
 141                iflag)))))
 142
 143 (in-package #:cl-user)
 144 #+#.(cl:if (cl:find-package '#:f2cl) '(and) '(or))
 145 (eval-when (:load-toplevel :compile-toplevel :execute)
 146   (setf (gethash 'fortran-to-lisp::dpcg fortran-to-lisp::*f2cl-function-info*)
 147           (fortran-to-lisp::make-f2cl-finfo
 148            :arg-types '((array fortran-to-lisp::integer4 (*)) (double-float)
 149                         (array double-float (*)) (array double-float (*))
 150                         (array double-float (*)) (array double-float (*))
 151                         (fortran-to-lisp::integer4) (fortran-to-lisp::integer4)
 152                         (double-float) (double-float)
 153                         (fortran-to-lisp::integer4) (fortran-to-lisp::integer4)
 154                         t t (fortran-to-lisp::integer4)
 155                         (array double-float (*)) (array double-float (*))
 156                         (array double-float (*)) (array double-float (*))
 157                         (fortran-to-lisp::integer4) (array double-float (*))
 158                         (array fortran-to-lisp::integer4 (*))
 159                         (array double-float (*)) (fortran-to-lisp::integer4))
 160            :return-values '(nil fortran-to-lisp::tn nil nil nil nil nil nil nil
 161                             fortran-to-lisp::hl0 nil nil nil nil
 162                             fortran-to-lisp::npsl nil nil nil nil
 163                             fortran-to-lisp::lpcg nil nil nil
 164                             fortran-to-lisp::iflag)
 165            :calls '(fortran-to-lisp::daxpy fortran-to-lisp::datp
 166                     fortran-to-lisp::ddot fortran-to-lisp::dcopy
 167                     fortran-to-lisp::dvnorm))))
 168