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(defun contrl 21 (xi xiold z dmz rhs delz deldmz dqz dqdmz g w v valstr slope scale 22 dscale accum ipvtg integs ipvtw nfxpnt fixpnt iflag fsub dfsub gsub 23 dgsub guess) 24 (declare (type (f2cl-lib:integer4) iflag nfxpnt) 25 (type (array f2cl-lib:integer4 (*)) ipvtw integs ipvtg) 26 (type (array double-float (*)) fixpnt accum dscale scale slope 27 valstr v w g dqdmz dqz deldmz delz 28 rhs dmz z xiold xi)) 29 (let ((colest-tolin 30 (make-array 40 31 :element-type 'double-float 32 :displaced-to (colest-part-0 *colest-common-block*) 33 :displaced-index-offset 120)) 34 (colest-ltol 35 (make-array 40 36 :element-type 'f2cl-lib:integer4 37 :displaced-to (colest-part-1 *colest-common-block*) 38 :displaced-index-offset 40))) 39 (symbol-macrolet ((precis (aref (colout-part-0 *colout-common-block*) 0)) 40 (iout (aref (colout-part-1 *colout-common-block*) 0)) 41 (iprint (aref (colout-part-1 *colout-common-block*) 1)) 42 (mstar (aref (colord-part-0 *colord-common-block*) 2)) 43 (kd (aref (colord-part-0 *colord-common-block*) 3)) 44 (n (aref (colapr-part-0 *colapr-common-block*) 0)) 45 (nold (aref (colapr-part-0 *colapr-common-block*) 1)) 46 (nmax (aref (colapr-part-0 *colapr-common-block*) 2)) 47 (nz (aref (colapr-part-0 *colapr-common-block*) 3)) 48 (ndmz (aref (colapr-part-0 *colapr-common-block*) 4)) 49 (mshnum (aref (colmsh-part-0 *colmsh-common-block*) 1)) 50 (mshlmt (aref (colmsh-part-0 *colmsh-common-block*) 2)) 51 (mshalt (aref (colmsh-part-0 *colmsh-common-block*) 3)) 52 (nonlin (aref (colnln-part-0 *colnln-common-block*) 0)) 53 (iter (aref (colnln-part-0 *colnln-common-block*) 1)) 54 (limit (aref (colnln-part-0 *colnln-common-block*) 2)) 55 (icare (aref (colnln-part-0 *colnln-common-block*) 3)) 56 (iguess (aref (colnln-part-0 *colnln-common-block*) 4)) 57 (tolin colest-tolin) 58 (ltol colest-ltol) 59 (ntol (aref (colest-part-1 *colest-common-block*) 80))) 60 (f2cl-lib:with-multi-array-data 61 ((xi double-float xi-%data% xi-%offset%) 62 (xiold double-float xiold-%data% xiold-%offset%) 63 (z double-float z-%data% z-%offset%) 64 (dmz double-float dmz-%data% dmz-%offset%) 65 (rhs double-float rhs-%data% rhs-%offset%) 66 (delz double-float delz-%data% delz-%offset%) 67 (deldmz double-float deldmz-%data% deldmz-%offset%) 68 (dqz double-float dqz-%data% dqz-%offset%) 69 (dqdmz double-float dqdmz-%data% dqdmz-%offset%) 70 (g double-float g-%data% g-%offset%) 71 (w double-float w-%data% w-%offset%) 72 (v double-float v-%data% v-%offset%) 73 (valstr double-float valstr-%data% valstr-%offset%) 74 (slope double-float slope-%data% slope-%offset%) 75 (scale double-float scale-%data% scale-%offset%) 76 (dscale double-float dscale-%data% dscale-%offset%) 77 (accum double-float accum-%data% accum-%offset%) 78 (fixpnt double-float fixpnt-%data% fixpnt-%offset%) 79 (ipvtg f2cl-lib:integer4 ipvtg-%data% ipvtg-%offset%) 80 (integs f2cl-lib:integer4 integs-%data% integs-%offset%) 81 (ipvtw f2cl-lib:integer4 ipvtw-%data% ipvtw-%offset%)) 82 (prog ((ifin 0) (lj 0) (j 0) (fact 0.0) (factor 0.0) (arg 0.0) 83 (anfix 0.0) (anorm 0.0) (ipred 0) (rlxold 0.0) (andif 0.0) 84 (anscl 0.0) (np1 0) (iz 0) (inz 0) (it 0) (ifrz 0) (rnold 0.0) 85 (ifreez 0) (relax 0.0) (rnorm 0.0) (msing 0) (noconv 0) (icor 0) 86 (iconv 0) (imesh 0) (i 0) (check 0.0) (lmtfrz 0) (rstart 0.0) 87 (relmin 0.0) (dummy (make-array 1 :element-type 'double-float))) 88 (declare (type (array double-float (1)) dummy) 89 (type double-float relmin rstart check rnorm relax rnold 90 anscl andif rlxold anorm anfix arg factor 91 fact) 92 (type (f2cl-lib:integer4) lmtfrz i imesh iconv icor noconv 93 msing ifreez ifrz it inz iz np1 94 ipred j lj ifin)) 95 (setf relmin 0.001) 96 (setf rstart 0.01) 97 (setf lmtfrz 4) 98 (setf check 0.0) 99 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) 100 ((> i ntol) nil) 101 (tagbody 102 label10 103 (setf check 104 (f2cl-lib:dmax1 (f2cl-lib:fref tolin (i) ((1 40))) 105 check)))) 106 (setf imesh 1) 107 (setf iconv 0) 108 (if (= nonlin 0) (setf iconv 1)) 109 (setf icor 0) 110 (setf noconv 0) 111 (setf msing 0) 112 label20 113 (setf iter 0) 114 (if (> nonlin 0) (go label50)) 115 (multiple-value-bind 116 (var-0 var-1 var-2 var-3 var-4 var-5 var-6 var-7 var-8 var-9 117 var-10 var-11 var-12 var-13 var-14 var-15 var-16 var-17 var-18 118 var-19 var-20 var-21) 119 (lsyslv msing xi xiold dummy dummy z dmz g w v rhs dummy integs 120 ipvtg ipvtw rnorm 0 fsub dfsub gsub dgsub guess) 121 (declare (ignore var-1 var-2 var-3 var-4 var-5 var-6 var-7 var-8 122 var-9 var-10 var-11 var-12 var-13 var-14 var-16 123 var-17 var-18 var-19 var-20 var-21)) 124 (setf msing var-0) 125 (setf rnorm var-15)) 126 (if (= msing 0) (go label400)) 127 label30 128 (if (< msing 0) (go label40)) 129 (if (< iprint 1) 130 (f2cl-lib:fformat iout 131 ("~%" "~%" 132 " A LOCAL ELIMINATION MATRIX IS SINGULAR " 133 "~%"))) 134 (go label460) 135 label40 136 (if (< iprint 1) 137 (f2cl-lib:fformat iout 138 ("~%" "~%" 139 " THE GLOBAL BVP-MATRIX IS SINGULAR " "~%"))) 140 (setf iflag 0) 141 (go end_label) 142 label50 143 (setf relax 1.0) 144 (if (or (= icare 1) (= icare -1)) (setf relax rstart)) 145 (if (= iconv 0) (go label160)) 146 (setf ifreez 0) 147 (multiple-value-bind 148 (var-0 var-1 var-2 var-3 var-4 var-5 var-6 var-7 var-8 var-9 149 var-10 var-11 var-12 var-13 var-14 var-15 var-16 var-17 var-18 150 var-19 var-20 var-21) 151 (lsyslv msing xi xiold z dmz delz deldmz g w v rhs dqdmz integs 152 ipvtg ipvtw rnold 1 fsub dfsub gsub dgsub guess) 153 (declare (ignore var-1 var-2 var-3 var-4 var-5 var-6 var-7 var-8 154 var-9 var-10 var-11 var-12 var-13 var-14 var-16 155 var-17 var-18 var-19 var-20 var-21)) 156 (setf msing var-0) 157 (setf rnold var-15)) 158 (if (< iprint 0) 159 (f2cl-lib:fformat iout 160 ("~%" " FIXED JACOBIAN ITERATIONS," "~%"))) 161 (if (< iprint 0) 162 (f2cl-lib:fformat iout 163 (" ITERATION = " 1 (("~3D")) " NORM (RHS) = " 164 1 (("~10,2,2,0,'*,,'DE")) "~%") 165 iter 166 rnold)) 167 (go label70) 168 label60 169 (if (< iprint 0) 170 (f2cl-lib:fformat iout 171 (" ITERATION = " 1 (("~3D")) " NORM (RHS) = " 172 1 (("~10,2,2,0,'*,,'DE")) "~%") 173 iter 174 rnorm)) 175 (setf rnold rnorm) 176 (multiple-value-bind 177 (var-0 var-1 var-2 var-3 var-4 var-5 var-6 var-7 var-8 var-9 178 var-10 var-11 var-12 var-13 var-14 var-15 var-16 var-17 var-18 179 var-19 var-20 var-21) 180 (lsyslv msing xi xiold z dmz delz deldmz g w v rhs dummy integs 181 ipvtg ipvtw rnorm (f2cl-lib:int-add 3 ifreez) fsub dfsub gsub 182 dgsub guess) 183 (declare (ignore var-1 var-2 var-3 var-4 var-5 var-6 var-7 var-8 184 var-9 var-10 var-11 var-12 var-13 var-14 var-16 185 var-17 var-18 var-19 var-20 var-21)) 186 (setf msing var-0) 187 (setf rnorm var-15)) 188 label70 189 (if (/= msing 0) (go label30)) 190 (if (= ifreez 1) (go label80)) 191 (setf iter (f2cl-lib:int-add iter 1)) 192 (setf ifrz 0) 193 label80 194 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) 195 ((> i nz) nil) 196 (tagbody 197 (setf (f2cl-lib:fref z-%data% (i) ((1 1)) z-%offset%) 198 (+ (f2cl-lib:fref z-%data% (i) ((1 1)) z-%offset%) 199 (f2cl-lib:fref delz-%data% 200 (i) 201 ((1 1)) 202 delz-%offset%))) 203 label90)) 204 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) 205 ((> i ndmz) nil) 206 (tagbody 207 (setf (f2cl-lib:fref dmz-%data% (i) ((1 1)) dmz-%offset%) 208 (+ (f2cl-lib:fref dmz-%data% (i) ((1 1)) dmz-%offset%) 209 (f2cl-lib:fref deldmz-%data% 210 (i) 211 ((1 1)) 212 deldmz-%offset%))) 213 label100)) 214 (multiple-value-bind 215 (var-0 var-1 var-2 var-3 var-4 var-5 var-6 var-7 var-8 var-9 216 var-10 var-11 var-12 var-13 var-14 var-15 var-16 var-17 var-18 217 var-19 var-20 var-21) 218 (lsyslv msing xi xiold z dmz delz deldmz g w v rhs dummy integs 219 ipvtg ipvtw rnorm 2 fsub dfsub gsub dgsub guess) 220 (declare (ignore var-1 var-2 var-3 var-4 var-5 var-6 var-7 var-8 221 var-9 var-10 var-11 var-12 var-13 var-14 var-16 222 var-17 var-18 var-19 var-20 var-21)) 223 (setf msing var-0) 224 (setf rnorm var-15)) 225 (if (< rnorm precis) (go label390)) 226 (if (> rnorm rnold) (go label130)) 227 (if (= ifreez 1) (go label110)) 228 (setf ifreez 1) 229 (go label60) 230 label110 231 (setf ifrz (f2cl-lib:int-add ifrz 1)) 232 (if (>= ifrz lmtfrz) (setf ifreez 0)) 233 (if (< rnold (* 4.0 rnorm)) (setf ifreez 0)) 234 (f2cl-lib:fdo (it 1 (f2cl-lib:int-add it 1)) 235 ((> it ntol) nil) 236 (tagbody 237 (setf inz (f2cl-lib:fref ltol (it) ((1 40)))) 238 (f2cl-lib:fdo (iz inz (f2cl-lib:int-add iz mstar)) 239 ((> iz nz) nil) 240 (tagbody 241 (if 242 (> 243 (f2cl-lib:dabs 244 (f2cl-lib:fref delz-%data% (iz) ((1 1)) delz-%offset%)) 245 (* (f2cl-lib:fref tolin (it) ((1 40))) 246 (+ 247 (f2cl-lib:dabs 248 (f2cl-lib:fref z-%data% (iz) ((1 1)) z-%offset%)) 249 1.0))) 250 (go label60)) 251 label120)))) 252 label120 253 (if (< iprint 1) 254 (f2cl-lib:fformat iout 255 ("~%" " CONVERGENCE AFTER" 1 (("~3D")) 256 " ITERATIONS" "~%" "~%") 257 iter)) 258 (go label400) 259 label130 260 (if (< iprint 0) 261 (f2cl-lib:fformat iout 262 (" ITERATION = " 1 (("~3D")) " NORM (RHS) = " 263 1 (("~10,2,2,0,'*,,'DE")) "~%") 264 iter 265 rnorm)) 266 (if (< iprint 0) 267 (f2cl-lib:fformat iout 268 ("~%" " SWITCH TO DAMPED NEWTON ITERATION," 269 "~%"))) 270 (setf iconv 0) 271 (setf relax rstart) 272 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) 273 ((> i nz) nil) 274 (tagbody 275 (setf (f2cl-lib:fref z-%data% (i) ((1 1)) z-%offset%) 276 (- (f2cl-lib:fref z-%data% (i) ((1 1)) z-%offset%) 277 (f2cl-lib:fref delz-%data% 278 (i) 279 ((1 1)) 280 delz-%offset%))) 281 label140)) 282 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) 283 ((> i ndmz) nil) 284 (tagbody 285 (setf (f2cl-lib:fref dmz-%data% (i) ((1 1)) dmz-%offset%) 286 (- (f2cl-lib:fref dmz-%data% (i) ((1 1)) dmz-%offset%) 287 (f2cl-lib:fref deldmz-%data% 288 (i) 289 ((1 1)) 290 deldmz-%offset%))) 291 label150)) 292 (setf np1 (f2cl-lib:int-add n 1)) 293 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) 294 ((> i np1) nil) 295 (tagbody 296 label155 297 (setf (f2cl-lib:fref xiold-%data% (i) ((1 1)) xiold-%offset%) 298 (f2cl-lib:fref xi-%data% (i) ((1 1)) xi-%offset%)))) 299 (setf nold n) 300 (setf iter 0) 301 label160 302 (if (< iprint 0) 303 (f2cl-lib:fformat iout 304 ("~%" " FULL DAMPED NEWTON ITERATION," "~%"))) 305 (multiple-value-bind 306 (var-0 var-1 var-2 var-3 var-4 var-5 var-6 var-7 var-8 var-9 307 var-10 var-11 var-12 var-13 var-14 var-15 var-16 var-17 var-18 308 var-19 var-20 var-21) 309 (lsyslv msing xi xiold z dmz delz deldmz g w v rhs dqdmz integs 310 ipvtg ipvtw rnold 1 fsub dfsub gsub dgsub guess) 311 (declare (ignore var-1 var-2 var-3 var-4 var-5 var-6 var-7 var-8 312 var-9 var-10 var-11 var-12 var-13 var-14 var-16 313 var-17 var-18 var-19 var-20 var-21)) 314 (setf msing var-0) 315 (setf rnold var-15)) 316 (if (/= msing 0) (go label30)) 317 (if (= iguess 1) (setf iguess 0)) 318 (skale n mstar kd z xi scale dscale) 319 (go label220) 320 label170 321 (setf rnold rnorm) 322 (if (>= iter limit) (go label430)) 323 (skale n mstar kd z xi scale dscale) 324 (setf anscl 0.0) 325 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) 326 ((> i nz) nil) 327 (tagbody 328 (setf anscl 329 (+ anscl 330 (expt 331 (* 332 (f2cl-lib:fref delz-%data% 333 (i) 334 ((1 1)) 335 delz-%offset%) 336 (f2cl-lib:fref scale-%data% 337 (i) 338 ((1 1)) 339 scale-%offset%)) 340 2))) 341 label180)) 342 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) 343 ((> i ndmz) nil) 344 (tagbody 345 (setf anscl 346 (+ anscl 347 (expt 348 (* 349 (f2cl-lib:fref deldmz-%data% 350 (i) 351 ((1 1)) 352 deldmz-%offset%) 353 (f2cl-lib:fref dscale-%data% 354 (i) 355 ((1 1)) 356 dscale-%offset%)) 357 2))) 358 label190)) 359 (setf anscl 360 (f2cl-lib:dsqrt 361 (/ anscl (f2cl-lib:dfloat (f2cl-lib:int-add nz ndmz))))) 362 (multiple-value-bind 363 (var-0 var-1 var-2 var-3 var-4 var-5 var-6 var-7 var-8 var-9 364 var-10 var-11 var-12 var-13 var-14 var-15 var-16 var-17 var-18 365 var-19 var-20 var-21) 366 (lsyslv msing xi xiold z dmz delz deldmz g w v rhs dummy integs 367 ipvtg ipvtw rnorm 3 fsub dfsub gsub dgsub guess) 368 (declare (ignore var-1 var-2 var-3 var-4 var-5 var-6 var-7 var-8 369 var-9 var-10 var-11 var-12 var-13 var-14 var-16 370 var-17 var-18 var-19 var-20 var-21)) 371 (setf msing var-0) 372 (setf rnorm var-15)) 373 (if (/= msing 0) (go label30)) 374 (setf andif 0.0) 375 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) 376 ((> i nz) nil) 377 (tagbody 378 (setf andif 379 (+ andif 380 (expt 381 (* 382 (- 383 (f2cl-lib:fref dqz-%data% (i) ((1 1)) dqz-%offset%) 384 (f2cl-lib:fref delz-%data% 385 (i) 386 ((1 1)) 387 delz-%offset%)) 388 (f2cl-lib:fref scale-%data% 389 (i) 390 ((1 1)) 391 scale-%offset%)) 392 2))) 393 label200)) 394 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) 395 ((> i ndmz) nil) 396 (tagbody 397 (setf andif 398 (+ andif 399 (expt 400 (* 401 (- 402 (f2cl-lib:fref dqdmz-%data% 403 (i) 404 ((1 1)) 405 dqdmz-%offset%) 406 (f2cl-lib:fref deldmz-%data% 407 (i) 408 ((1 1)) 409 deldmz-%offset%)) 410 (f2cl-lib:fref dscale-%data% 411 (i) 412 ((1 1)) 413 dscale-%offset%)) 414 2))) 415 label210)) 416 (setf andif 417 (f2cl-lib:dsqrt 418 (+ (/ andif (f2cl-lib:dfloat (f2cl-lib:int-add nz ndmz))) 419 precis))) 420 (setf relax (/ (* relax anscl) andif)) 421 (if (> relax 1.0) (setf relax 1.0)) 422 label220 423 (setf rlxold relax) 424 (setf ipred 1) 425 (setf iter (f2cl-lib:int-add iter 1)) 426 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) 427 ((> i nz) nil) 428 (tagbody 429 (setf (f2cl-lib:fref z-%data% (i) ((1 1)) z-%offset%) 430 (+ (f2cl-lib:fref z-%data% (i) ((1 1)) z-%offset%) 431 (* relax 432 (f2cl-lib:fref delz-%data% 433 (i) 434 ((1 1)) 435 delz-%offset%)))) 436 label230)) 437 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) 438 ((> i ndmz) nil) 439 (tagbody 440 (setf (f2cl-lib:fref dmz-%data% (i) ((1 1)) dmz-%offset%) 441 (+ (f2cl-lib:fref dmz-%data% (i) ((1 1)) dmz-%offset%) 442 (* relax 443 (f2cl-lib:fref deldmz-%data% 444 (i) 445 ((1 1)) 446 deldmz-%offset%)))) 447 label240)) 448 label250 449 (multiple-value-bind 450 (var-0 var-1 var-2 var-3 var-4 var-5 var-6 var-7 var-8 var-9 451 var-10 var-11 var-12 var-13 var-14 var-15 var-16 var-17 var-18 452 var-19 var-20 var-21) 453 (lsyslv msing xi xiold z dmz dqz dqdmz g w v rhs dummy integs 454 ipvtg ipvtw rnorm 2 fsub dfsub gsub dgsub guess) 455 (declare (ignore var-1 var-2 var-3 var-4 var-5 var-6 var-7 var-8 456 var-9 var-10 var-11 var-12 var-13 var-14 var-16 457 var-17 var-18 var-19 var-20 var-21)) 458 (setf msing var-0) 459 (setf rnorm var-15)) 460 (multiple-value-bind 461 (var-0 var-1 var-2 var-3 var-4 var-5 var-6 var-7 var-8 var-9 462 var-10 var-11 var-12 var-13 var-14 var-15 var-16 var-17 var-18 463 var-19 var-20 var-21) 464 (lsyslv msing xi xiold z dmz dqz dqdmz g w v rhs dummy integs 465 ipvtg ipvtw rnorm 4 fsub dfsub gsub dgsub guess) 466 (declare (ignore var-1 var-2 var-3 var-4 var-5 var-6 var-7 var-8 467 var-9 var-10 var-11 var-12 var-13 var-14 var-16 468 var-17 var-18 var-19 var-20 var-21)) 469 (setf msing var-0) 470 (setf rnorm var-15)) 471 (setf anorm 0.0) 472 (setf anfix 0.0) 473 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) 474 ((> i nz) nil) 475 (tagbody 476 (setf anorm 477 (+ anorm 478 (expt 479 (* 480 (f2cl-lib:fref delz-%data% 481 (i) 482 ((1 1)) 483 delz-%offset%) 484 (f2cl-lib:fref scale-%data% 485 (i) 486 ((1 1)) 487 scale-%offset%)) 488 2))) 489 (setf anfix 490 (+ anfix 491 (expt 492 (* 493 (f2cl-lib:fref dqz-%data% (i) ((1 1)) dqz-%offset%) 494 (f2cl-lib:fref scale-%data% 495 (i) 496 ((1 1)) 497 scale-%offset%)) 498 2))) 499 label260)) 500 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) 501 ((> i ndmz) nil) 502 (tagbody 503 (setf anorm 504 (+ anorm 505 (expt 506 (* 507 (f2cl-lib:fref deldmz-%data% 508 (i) 509 ((1 1)) 510 deldmz-%offset%) 511 (f2cl-lib:fref dscale-%data% 512 (i) 513 ((1 1)) 514 dscale-%offset%)) 515 2))) 516 (setf anfix 517 (+ anfix 518 (expt 519 (* 520 (f2cl-lib:fref dqdmz-%data% 521 (i) 522 ((1 1)) 523 dqdmz-%offset%) 524 (f2cl-lib:fref dscale-%data% 525 (i) 526 ((1 1)) 527 dscale-%offset%)) 528 2))) 529 label270)) 530 (setf anorm 531 (f2cl-lib:dsqrt 532 (/ anorm (f2cl-lib:dfloat (f2cl-lib:int-add nz ndmz))))) 533 (setf anfix 534 (f2cl-lib:dsqrt 535 (/ anfix (f2cl-lib:dfloat (f2cl-lib:int-add nz ndmz))))) 536 (if (= icor 1) (go label280)) 537 (if (< iprint 0) 538 (f2cl-lib:fformat iout 539 (" ITERATION = " 1 (("~3D")) 540 " RELAXATION FACTOR = " 1 541 (("~10,2,2,0,'*,,'DE")) "~%" 542 " NORM OF SCALED RHS CHANGES FROM " 1 543 (("~10,2,2,0,'*,,'DE")) " TO" 1 544 (("~10,2,2,0,'*,,'DE")) "~%" 545 " NORM OF RHS CHANGES FROM " 1 546 (("~10,2,2,0,'*,,'DE")) " TO" 1 547 (("~10,2,2,0,'*,,'DE")) 1 548 (("~10,2,2,0,'*,,'DE")) "~%") 549 iter 550 relax 551 anorm 552 anfix 553 rnold 554 rnorm)) 555 (go label290) 556 label280 557 (if (< iprint 0) 558 (f2cl-lib:fformat iout 559 (" RELAXATION FACTOR CORRECTED TO RELAX = " 1 560 (("~10,2,2,0,'*,,'DE")) "~%" 561 " NORM OF SCALED RHS CHANGES FROM " 1 562 (("~10,2,2,0,'*,,'DE")) " TO" 1 563 (("~10,2,2,0,'*,,'DE")) "~%" 564 " NORM OF RHS CHANGES FROM " 1 565 (("~10,2,2,0,'*,,'DE")) " TO" 1 566 (("~10,2,2,0,'*,,'DE")) 1 567 (("~10,2,2,0,'*,,'DE")) "~%") 568 relax 569 anorm 570 anfix 571 rnold 572 rnorm)) 573 label290 574 (setf icor 0) 575 (if (or (< anfix precis) (< rnorm precis)) (go label390)) 576 (if (> anfix anorm) (go label300)) 577 (if (<= anfix check) (go label350)) 578 (if (/= ipred 1) (go label170)) 579 label300 580 (if (>= iter limit) (go label430)) 581 (setf ipred 0) 582 (setf arg (+ (/ (- (/ anfix anorm) 1.0) relax) 1.0)) 583 (if (< arg 0.0) (go label170)) 584 (if (<= arg (+ (* 0.25 relax) (* 0.125 (expt relax 2)))) 585 (go label310)) 586 (setf factor (- (f2cl-lib:dsqrt (+ 1.0 (* 8.0 arg))) 1.0)) 587 (if (< (f2cl-lib:dabs (- factor 1.0)) (* 0.1 factor)) (go label170)) 588 (if (< factor 0.5) (setf factor 0.5)) 589 (setf relax (/ relax factor)) 590 (go label320) 591 label310 592 (if (>= relax 0.9) (go label170)) 593 (setf relax 1.0) 594 label320 595 (setf icor 1) 596 (if (< relax relmin) (go label440)) 597 (setf fact (- relax rlxold)) 598 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) 599 ((> i nz) nil) 600 (tagbody 601 (setf (f2cl-lib:fref z-%data% (i) ((1 1)) z-%offset%) 602 (+ (f2cl-lib:fref z-%data% (i) ((1 1)) z-%offset%) 603 (* fact 604 (f2cl-lib:fref delz-%data% 605 (i) 606 ((1 1)) 607 delz-%offset%)))) 608 label330)) 609 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) 610 ((> i ndmz) nil) 611 (tagbody 612 (setf (f2cl-lib:fref dmz-%data% (i) ((1 1)) dmz-%offset%) 613 (+ (f2cl-lib:fref dmz-%data% (i) ((1 1)) dmz-%offset%) 614 (* fact 615 (f2cl-lib:fref deldmz-%data% 616 (i) 617 ((1 1)) 618 deldmz-%offset%)))) 619 label340)) 620 (setf rlxold relax) 621 (go label250) 622 label350 623 (f2cl-lib:fdo (it 1 (f2cl-lib:int-add it 1)) 624 ((> it ntol) nil) 625 (tagbody 626 (setf inz (f2cl-lib:fref ltol (it) ((1 40)))) 627 (f2cl-lib:fdo (iz inz (f2cl-lib:int-add iz mstar)) 628 ((> iz nz) nil) 629 (tagbody 630 (if 631 (> 632 (f2cl-lib:dabs 633 (f2cl-lib:fref dqz-%data% (iz) ((1 1)) dqz-%offset%)) 634 (* (f2cl-lib:fref tolin (it) ((1 40))) 635 (+ 636 (f2cl-lib:dabs 637 (f2cl-lib:fref z-%data% (iz) ((1 1)) z-%offset%)) 638 1.0))) 639 (go label170)) 640 label360)))) 641 label360 642 (if (< iprint 1) 643 (f2cl-lib:fformat iout 644 ("~%" " CONVERGENCE AFTER" 1 (("~3D")) 645 " ITERATIONS" "~%" "~%") 646 iter)) 647 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) 648 ((> i nz) nil) 649 (tagbody 650 (setf (f2cl-lib:fref z-%data% (i) ((1 1)) z-%offset%) 651 (+ (f2cl-lib:fref z-%data% (i) ((1 1)) z-%offset%) 652 (f2cl-lib:fref dqz-%data% (i) ((1 1)) dqz-%offset%))) 653 label370)) 654 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) 655 ((> i ndmz) nil) 656 (tagbody 657 (setf (f2cl-lib:fref dmz-%data% (i) ((1 1)) dmz-%offset%) 658 (+ (f2cl-lib:fref dmz-%data% (i) ((1 1)) dmz-%offset%) 659 (f2cl-lib:fref dqdmz-%data% 660 (i) 661 ((1 1)) 662 dqdmz-%offset%))) 663 label380)) 664 label390 665 (if (and (or (< anfix precis) (< rnorm precis)) (< iprint 1)) 666 (f2cl-lib:fformat iout 667 ("~%" " CONVERGENCE AFTER" 1 (("~3D")) 668 " ITERATIONS" "~%" "~%") 669 iter)) 670 (setf iconv 1) 671 (if (= icare -1) (setf icare 0)) 672 label400 673 (if (>= iprint 0) (go label420)) 674 (f2cl-lib:fdo (j 1 (f2cl-lib:int-add j 1)) 675 ((> j mstar) nil) 676 (tagbody 677 (f2cl-lib:fformat iout 678 (" MESH VALUES FOR Z(" 1 (("~2D")) ")," "~%") 679 j) 680 label410 681 (f2cl-lib:fformat iout 682 (" " 8 (("~15,7,2,0,'*,,'DE")) "~%") 683 (do ((lj j (f2cl-lib:int-add lj mstar)) 684 (%ret nil)) 685 ((> lj nz) (nreverse %ret)) 686 (declare (type f2cl-lib:integer4 lj)) 687 (push 688 (f2cl-lib:fref z-%data% 689 (lj) 690 ((1 1)) 691 z-%offset%) 692 %ret))))) 693 label420 694 (setf ifin 1) 695 (if (= imesh 2) 696 (multiple-value-bind (var-0 var-1 var-2 var-3 var-4) 697 (errchk xi z dmz valstr ifin) 698 (declare (ignore var-0 var-1 var-2 var-3)) 699 (setf ifin var-4))) 700 (if (or (= imesh 1) (and (= ifin 0) (/= icare 2))) (go label460)) 701 (setf iflag 1) 702 (go end_label) 703 label430 704 (if (< iprint 1) 705 (f2cl-lib:fformat iout 706 ("~%" " NO CONVERGENCE AFTER " 1 (("~3D")) 707 " ITERATIONS" "~%" "~%") 708 iter)) 709 (go label450) 710 label440 711 (if (< iprint 1) 712 (f2cl-lib:fformat iout 713 ("~%" " NO CONVERGENCE. RELAXATION FACTOR =" 1 714 (("~10,3,2,0,'*,,'DE")) 715 " IS TOO SMALL (LESS THAN" 1 716 (("~10,3,2,0,'*,,'DE")) ")" "~%" "~%") 717 relax 718 relmin)) 719 label450 720 (setf iflag -2) 721 (setf noconv (f2cl-lib:int-add noconv 1)) 722 (if (and (= icare 2) (> noconv 1)) (go end_label)) 723 (if (= icare 0) (setf icare -1)) 724 label460 725 (setf np1 (f2cl-lib:int-add n 1)) 726 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) 727 ((> i np1) nil) 728 (tagbody 729 label470 730 (setf (f2cl-lib:fref xiold-%data% (i) ((1 1)) xiold-%offset%) 731 (f2cl-lib:fref xi-%data% (i) ((1 1)) xi-%offset%)))) 732 (setf nold n) 733 (setf imesh 1) 734 (if (or (= iconv 0) (>= mshnum mshlmt) (>= mshalt mshlmt)) 735 (setf imesh 2)) 736 (if (and (>= mshalt mshlmt) (< mshnum mshlmt)) (setf mshalt 1)) 737 (multiple-value-bind 738 (var-0 var-1 var-2 var-3 var-4 var-5 var-6 var-7 var-8 var-9) 739 (newmsh imesh xi xiold z dmz valstr slope accum nfxpnt fixpnt) 740 (declare (ignore var-1 var-2 var-3 var-4 var-5 var-6 var-7 var-8 741 var-9)) 742 (setf imesh var-0)) 743 (if (<= n nmax) (go label480)) 744 (setf n (the f2cl-lib:integer4 (truncate n 2))) 745 (setf iflag -1) 746 (if (and (= iconv 0) (< iprint 1)) 747 (f2cl-lib:fformat iout (" (NO CONVERGENCE)" "~%"))) 748 (if (and (= iconv 1) (< iprint 1)) 749 (f2cl-lib:fformat iout 750 (" (PROBABLY TOLERANCES TOO STRINGENT, OR NMAX TOO " 751 "SMALL)" "~%"))) 752 (go end_label) 753 label480 754 (if (= iconv 0) (setf imesh 1)) 755 (if (= icare 1) (setf iconv 0)) 756 (go label20) 757 end_label 758 (return 759 (values nil 760 nil 761 nil 762 nil 763 nil 764 nil 765 nil 766 nil 767 nil 768 nil 769 nil 770 nil 771 nil 772 nil 773 nil 774 nil 775 nil 776 nil 777 nil 778 nil 779 nil 780 nil 781 iflag 782 nil 783 nil 784 nil 785 nil 786 nil))))))) 787 788(in-package #-gcl #:cl-user #+gcl "CL-USER") 789#+#.(cl:if (cl:find-package '#:f2cl) '(and) '(or)) 790(eval-when (:load-toplevel :compile-toplevel :execute) 791 (setf (gethash 'fortran-to-lisp::contrl 792 fortran-to-lisp::*f2cl-function-info*) 793 (fortran-to-lisp::make-f2cl-finfo 794 :arg-types '((array double-float (1)) (array double-float (1)) 795 (array double-float (1)) (array double-float (1)) 796 (array double-float (1)) (array double-float (1)) 797 (array double-float (1)) (array double-float (1)) 798 (array double-float (1)) (array double-float (1)) 799 (array double-float (1)) (array double-float (1)) 800 (array double-float (1)) (array double-float (1)) 801 (array double-float (1)) (array double-float (1)) 802 (array double-float (1)) 803 (array fortran-to-lisp::integer4 (1)) 804 (array fortran-to-lisp::integer4 (1)) 805 (array fortran-to-lisp::integer4 (1)) 806 (fortran-to-lisp::integer4) (array double-float (1)) 807 (fortran-to-lisp::integer4) t t t t t) 808 :return-values '(nil nil nil nil nil nil nil nil nil nil nil nil nil 809 nil nil nil nil nil nil nil nil nil 810 fortran-to-lisp::iflag nil nil nil nil nil) 811 :calls '(fortran-to-lisp::newmsh fortran-to-lisp::errchk 812 fortran-to-lisp::skale fortran-to-lisp::lsyslv)))) 813 814