share/colnew/lisp/consts.lisp

   1 ;;; Compiled by f2cl version:
   2 ;;; ("f2cl1.l,v 1.221 2010/05/26 19:25:52 rtoy Exp $"
   3 ;;;  "f2cl2.l,v 1.37 2008/02/22 22:19:33 rtoy Exp $"
   4 ;;;  "f2cl3.l,v 1.6 2008/02/22 22:19:33 rtoy Exp $"
   5 ;;;  "f2cl4.l,v 1.7 2008/02/22 22:19:34 rtoy Exp $"
   6 ;;;  "f2cl5.l,v 1.204 2010/02/23 05:21:30 rtoy Exp $"
   7 ;;;  "f2cl6.l,v 1.48 2008/08/24 00:56:27 rtoy Exp $"
   8 ;;;  "macros.l,v 1.114 2010/05/17 01:42:14 rtoy Exp $")
   9
  10 ;;; Using Lisp CMU Common Lisp CVS Head 2010-05-25 18:21:07 (20A 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 double-float))
  16
  17 (in-package :colnew)
  18
  19
  20 (let ((cnsts1
  21        (make-array 28
  22                    :element-type 'double-float
  23                    :initial-contents '(0.25 0.0625 0.072169 0.018342 0.019065
  24                                        0.05819 0.0054658 0.005337 0.01889
  25                                        0.027792 0.0016095 0.0014964 0.0075938
  26                                        0.0057573 0.018342 0.004673 4.15e-4
  27                                        0.001919 0.001468 0.006371 0.00461
  28                                        1.342e-4 1.138e-4 4.889e-4 4.177e-4
  29                                        0.001374 0.001654 0.002863)))
  30       (cnsts2
  31        (make-array 28
  32                    :element-type 'double-float
  33                    :initial-contents '(0.125 0.002604 0.008019 2.17e-5 7.453e-5
  34                                        5.208e-4 9.689e-8 3.689e-7 3.1e-6
  35                                        2.451e-5 2.691e-10 1.12e-9 1.076e-8
  36                                        9.405e-8 1.033e-6 5.097e-13 2.29e-12
  37                                        2.446e-11 2.331e-10 2.936e-9 3.593e-8
  38                                        7.001e-16 3.363e-15 3.921e-14 4.028e-13
  39                                        5.646e-12 7.531e-11 1.129e-9))))
  40   (declare (type (array double-float (28)) cnsts1 cnsts2))
  41   (defun consts (k rho coef)
  42     (declare (type (array double-float (*)) coef)
  43              (type (array double-float (*)) rho)
  44              (type (f2cl-lib:integer4) k))
  45     (let ((colord-m
  46            (make-array 20
  47                        :element-type 'f2cl-lib:integer4
  48                        :displaced-to (colord-part-0 *colord-common-block*)
  49                        :displaced-index-offset 5))
  50           (colbas-b
  51            (make-array 28
  52                        :element-type 'double-float
  53                        :displaced-to (colbas-part-0 *colbas-common-block*)
  54                        :displaced-index-offset 0))
  55           (colbas-acol
  56            (make-array 196
  57                        :element-type 'double-float
  58                        :displaced-to (colbas-part-0 *colbas-common-block*)
  59                        :displaced-index-offset 28))
  60           (colbas-asave
  61            (make-array 112
  62                        :element-type 'double-float
  63                        :displaced-to (colbas-part-0 *colbas-common-block*)
  64                        :displaced-index-offset 224))
  65           (colest-wgtmsh
  66            (make-array 40
  67                        :element-type 'double-float
  68                        :displaced-to (colest-part-0 *colest-common-block*)
  69                        :displaced-index-offset 40))
  70           (colest-wgterr
  71            (make-array 40
  72                        :element-type 'double-float
  73                        :displaced-to (colest-part-0 *colest-common-block*)
  74                        :displaced-index-offset 80))
  75           (colest-tolin
  76            (make-array 40
  77                        :element-type 'double-float
  78                        :displaced-to (colest-part-0 *colest-common-block*)
  79                        :displaced-index-offset 120))
  80           (colest-root
  81            (make-array 40
  82                        :element-type 'double-float
  83                        :displaced-to (colest-part-0 *colest-common-block*)
  84                        :displaced-index-offset 160))
  85           (colest-jtol
  86            (make-array 40
  87                        :element-type 'f2cl-lib:integer4
  88                        :displaced-to (colest-part-1 *colest-common-block*)
  89                        :displaced-index-offset 0))
  90           (colest-ltol
  91            (make-array 40
  92                        :element-type 'f2cl-lib:integer4
  93                        :displaced-to (colest-part-1 *colest-common-block*)
  94                        :displaced-index-offset 40)))
  95       (symbol-macrolet ((ncomp (aref (colord-part-0 *colord-common-block*) 1))
  96                         (mmax (aref (colord-part-0 *colord-common-block*) 4))
  97                         (m colord-m)
  98                         (b colbas-b)
  99                         (acol colbas-acol)
 100                         (asave colbas-asave)
 101                         (wgtmsh colest-wgtmsh)
 102                         (wgterr colest-wgterr)
 103                         (tolin colest-tolin)
 104                         (root colest-root)
 105                         (jtol colest-jtol)
 106                         (ltol colest-ltol)
 107                         (ntol (aref (colest-part-1 *colest-common-block*) 80)))
 108         (f2cl-lib:with-multi-array-data
 109             ((rho double-float rho-%data% rho-%offset%)
 110              (coef double-float coef-%data% coef-%offset%))
 111           (prog ((ltoli 0) (i 0) (mtot 0) (jcomp 0) (l 0) (mj 0) (j 0) (iz 0)
 112                  (koff 0) (dummy (make-array 1 :element-type 'double-float)))
 113             (declare (type (array double-float (1)) dummy)
 114                      (type (f2cl-lib:integer4) koff iz j mj l jcomp mtot i
 115                                                ltoli))
 116             (setf koff (the f2cl-lib:integer4 (truncate (* k (+ k 1)) 2)))
 117             (setf iz 1)
 118             (f2cl-lib:fdo (j 1 (f2cl-lib:int-add j 1))
 119                           ((> j ncomp) nil)
 120               (tagbody
 121                 (setf mj (f2cl-lib:fref m (j) ((1 20))))
 122                 (f2cl-lib:fdo (l 1 (f2cl-lib:int-add l 1))
 123                               ((> l mj) nil)
 124                   (tagbody
 125                     (setf (f2cl-lib:fref wgterr (iz) ((1 40)))
 126                             (f2cl-lib:fref cnsts1
 127                                            ((f2cl-lib:int-add
 128                                              (f2cl-lib:int-sub koff mj)
 129                                              l))
 130                                            ((1 28))))
 131                     (setf iz (f2cl-lib:int-add iz 1))
 132                    label10))))
 133            label10
 134             (setf jcomp 1)
 135             (setf mtot (f2cl-lib:fref m (1) ((1 20))))
 136             (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1))
 137                           ((> i ntol) nil)
 138               (tagbody
 139                 (setf ltoli (f2cl-lib:fref ltol (i) ((1 40))))
 140                label20
 141                 (if (<= ltoli mtot) (go label30))
 142                 (setf jcomp (f2cl-lib:int-add jcomp 1))
 143                 (setf mtot
 144                         (f2cl-lib:int-add mtot
 145                                           (f2cl-lib:fref m (jcomp) ((1 20)))))
 146                 (go label20)
 147                label30
 148                 (setf (f2cl-lib:fref jtol (i) ((1 40))) jcomp)
 149                 (setf (f2cl-lib:fref wgtmsh (i) ((1 40)))
 150                         (/
 151                          (* 10.0
 152                             (f2cl-lib:fref cnsts2
 153                                            ((f2cl-lib:int-sub
 154                                              (f2cl-lib:int-add koff ltoli)
 155                                              mtot))
 156                                            ((1 28))))
 157                          (f2cl-lib:fref tolin (i) ((1 40)))))
 158                 (setf (f2cl-lib:fref root (i) ((1 40)))
 159                         (/ 1.0
 160                            (f2cl-lib:dfloat
 161                             (f2cl-lib:int-add
 162                              (f2cl-lib:int-sub (f2cl-lib:int-add k mtot) ltoli)
 163                              1))))
 164                label40))
 165             (f2cl-lib:computed-goto
 166              (label50 label60 label70 label80 label90 label100 label110)
 167              k)
 168            label50
 169             (setf (f2cl-lib:fref rho-%data% (1) ((1 7)) rho-%offset%) 0.0)
 170             (go label120)
 171            label60
 172             (setf (f2cl-lib:fref rho-%data% (2) ((1 7)) rho-%offset%)
 173                     0.5773502691896257)
 174             (setf (f2cl-lib:fref rho-%data% (1) ((1 7)) rho-%offset%)
 175                     (- (f2cl-lib:fref rho-%data% (2) ((1 7)) rho-%offset%)))
 176             (go label120)
 177            label70
 178             (setf (f2cl-lib:fref rho-%data% (3) ((1 7)) rho-%offset%)
 179                     0.7745966692414834)
 180             (setf (f2cl-lib:fref rho-%data% (2) ((1 7)) rho-%offset%) 0.0)
 181             (setf (f2cl-lib:fref rho-%data% (1) ((1 7)) rho-%offset%)
 182                     (- (f2cl-lib:fref rho-%data% (3) ((1 7)) rho-%offset%)))
 183             (go label120)
 184            label80
 185             (setf (f2cl-lib:fref rho-%data% (4) ((1 7)) rho-%offset%)
 186                     0.8611363115940526)
 187             (setf (f2cl-lib:fref rho-%data% (3) ((1 7)) rho-%offset%)
 188                     0.33998104358485626)
 189             (setf (f2cl-lib:fref rho-%data% (2) ((1 7)) rho-%offset%)
 190                     (- (f2cl-lib:fref rho-%data% (3) ((1 7)) rho-%offset%)))
 191             (setf (f2cl-lib:fref rho-%data% (1) ((1 7)) rho-%offset%)
 192                     (- (f2cl-lib:fref rho-%data% (4) ((1 7)) rho-%offset%)))
 193             (go label120)
 194            label90
 195             (setf (f2cl-lib:fref rho-%data% (5) ((1 7)) rho-%offset%)
 196                     0.906179845938664)
 197             (setf (f2cl-lib:fref rho-%data% (4) ((1 7)) rho-%offset%)
 198                     0.5384693101056831)
 199             (setf (f2cl-lib:fref rho-%data% (3) ((1 7)) rho-%offset%) 0.0)
 200             (setf (f2cl-lib:fref rho-%data% (2) ((1 7)) rho-%offset%)
 201                     (- (f2cl-lib:fref rho-%data% (4) ((1 7)) rho-%offset%)))
 202             (setf (f2cl-lib:fref rho-%data% (1) ((1 7)) rho-%offset%)
 203                     (- (f2cl-lib:fref rho-%data% (5) ((1 7)) rho-%offset%)))
 204             (go label120)
 205            label100
 206             (setf (f2cl-lib:fref rho-%data% (6) ((1 7)) rho-%offset%)
 207                     0.932469514203152)
 208             (setf (f2cl-lib:fref rho-%data% (5) ((1 7)) rho-%offset%)
 209                     0.6612093864662645)
 210             (setf (f2cl-lib:fref rho-%data% (4) ((1 7)) rho-%offset%)
 211                     0.2386191860831969)
 212             (setf (f2cl-lib:fref rho-%data% (3) ((1 7)) rho-%offset%)
 213                     (- (f2cl-lib:fref rho-%data% (4) ((1 7)) rho-%offset%)))
 214             (setf (f2cl-lib:fref rho-%data% (2) ((1 7)) rho-%offset%)
 215                     (- (f2cl-lib:fref rho-%data% (5) ((1 7)) rho-%offset%)))
 216             (setf (f2cl-lib:fref rho-%data% (1) ((1 7)) rho-%offset%)
 217                     (- (f2cl-lib:fref rho-%data% (6) ((1 7)) rho-%offset%)))
 218             (go label120)
 219            label110
 220             (setf (f2cl-lib:fref rho-%data% (7) ((1 7)) rho-%offset%)
 221                     0.9491079912342758)
 222             (setf (f2cl-lib:fref rho-%data% (6) ((1 7)) rho-%offset%)
 223                     0.7415311855993945)
 224             (setf (f2cl-lib:fref rho-%data% (5) ((1 7)) rho-%offset%)
 225                     0.4058451513773972)
 226             (setf (f2cl-lib:fref rho-%data% (4) ((1 7)) rho-%offset%) 0.0)
 227             (setf (f2cl-lib:fref rho-%data% (3) ((1 7)) rho-%offset%)
 228                     (- (f2cl-lib:fref rho-%data% (5) ((1 7)) rho-%offset%)))
 229             (setf (f2cl-lib:fref rho-%data% (2) ((1 7)) rho-%offset%)
 230                     (- (f2cl-lib:fref rho-%data% (6) ((1 7)) rho-%offset%)))
 231             (setf (f2cl-lib:fref rho-%data% (1) ((1 7)) rho-%offset%)
 232                     (- (f2cl-lib:fref rho-%data% (7) ((1 7)) rho-%offset%)))
 233            label120
 234             (f2cl-lib:fdo (j 1 (f2cl-lib:int-add j 1))
 235                           ((> j k) nil)
 236               (tagbody
 237                 (setf (f2cl-lib:fref rho-%data% (j) ((1 7)) rho-%offset%)
 238                         (* 0.5
 239                            (+ 1.0
 240                               (f2cl-lib:fref rho-%data%
 241                                              (j)
 242                                              ((1 7))
 243                                              rho-%offset%))))
 244                label130))
 245             (f2cl-lib:fdo (j 1 (f2cl-lib:int-add j 1))
 246                           ((> j k) nil)
 247               (tagbody
 248                 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1))
 249                               ((> i k) nil)
 250                   (tagbody
 251                    label135
 252                     (setf (f2cl-lib:fref coef-%data%
 253                                          (i j)
 254                                          ((1 k) (1 1))
 255                                          coef-%offset%)
 256                             0.0)))
 257                 (setf (f2cl-lib:fref coef-%data%
 258                                      (j j)
 259                                      ((1 k) (1 1))
 260                                      coef-%offset%)
 261                         1.0)
 262                 (vmonde rho
 263                  (f2cl-lib:array-slice coef double-float (1 j) ((1 k) (1 1)))
 264                  k)
 265                label140))
 266             (rkbas 1.0 coef k mmax b dummy 0)
 267             (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1))
 268                           ((> i k) nil)
 269               (tagbody
 270                 (rkbas (f2cl-lib:fref rho-%data% (i) ((1 7)) rho-%offset%) coef
 271                  k mmax
 272                  (f2cl-lib:array-slice acol double-float (1 i) ((1 28) (1 7)))
 273                  dummy 0)
 274                label150))
 275             (rkbas (/ 1.0 6.0) coef k mmax
 276              (f2cl-lib:array-slice asave double-float (1 1) ((1 28) (1 4)))
 277              dummy 0)
 278             (rkbas (/ 1.0 3.0) coef k mmax
 279              (f2cl-lib:array-slice asave double-float (1 2) ((1 28) (1 4)))
 280              dummy 0)
 281             (rkbas (/ 2.0 3.0) coef k mmax
 282              (f2cl-lib:array-slice asave double-float (1 3) ((1 28) (1 4)))
 283              dummy 0)
 284             (rkbas (/ 5.0 6.0) coef k mmax
 285              (f2cl-lib:array-slice asave double-float (1 4) ((1 28) (1 4)))
 286              dummy 0)
 287             (go end_label)
 288            end_label
 289             (return (values nil nil nil))))))))
 290
 291 (in-package #-gcl #:cl-user #+gcl "CL-USER")
 292 #+#.(cl:if (cl:find-package '#:f2cl) '(and) '(or))
 293 (eval-when (:load-toplevel :compile-toplevel :execute)
 294   (setf (gethash 'fortran-to-lisp::consts
 295                  fortran-to-lisp::*f2cl-function-info*)
 296           (fortran-to-lisp::make-f2cl-finfo
 297            :arg-types '((fortran-to-lisp::integer4) (array double-float (7))
 298                         (array double-float (*)))
 299            :return-values '(nil nil nil)
 300            :calls '(fortran-to-lisp::rkbas fortran-to-lisp::vmonde))))
 301