?  search:
summary | log | graphiclog1 | graphiclog2 | commit | commitdiff | tree | refs | edit | fork
blame | history | raw | HEAD
Remove the obsolete DEFMTRFUN-EXTERNAL macro
[maxima.git] / share / hompack / lisp / pcgns.lisp
blob0cfd1466c066e8dd6e6dbf8fabd955af12cdf870
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-2020-04 (21D 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 "HOMPACK")
18
19
20 (defun pcgns (nn aa lenaa maxa pp rho start work iflag)
21 (declare (type (array f2cl-lib:integer4 (*)) maxa)
22 (type (array double-float (*)) work start rho pp aa)
23 (type (f2cl-lib:integer4) iflag lenaa nn))
24 (f2cl-lib:with-multi-array-data
25 ((aa double-float aa-%data% aa-%offset%)
26 (pp double-float pp-%data% pp-%offset%)
27 (rho double-float rho-%data% rho-%offset%)
28 (start double-float start-%data% start-%offset%)
29 (work double-float work-%data% work-%offset%)
30 (maxa f2cl-lib:integer4 maxa-%data% maxa-%offset%))
31 (prog ((stillu nil) (stillb nil) (ab 0.0) (au 0.0) (bb 0.0) (bu 0.0)
32 (dznrm 0.0) (pbnprd 0.0) (punprd 0.0) (rbnprd 0.0) (rbtol 0.0)
33 (rnprd 0.0) (runprd 0.0) (rutol 0.0) (startk 0.0) (temp 0.0)
34 (unrm 0.0) (zlen 0.0) (ztol 0.0) (imax 0) (ind 0) (j 0) (k 0)
35 (np1 0) (np2 0) (n2p3 0) (n3p4 0) (n4p5 0) (n5p6 0))
36 (declare (type (f2cl-lib:integer4) n5p6 n4p5 n3p4 n2p3 np2 np1 k j ind
37 imax)
38 (type (double-float) ztol zlen unrm temp startk rutol runprd
39 rnprd rbtol rbnprd punprd pbnprd dznrm bu
40 bb au ab)
41 (type f2cl-lib:logical stillb stillu))
42 (setf np1 (f2cl-lib:int-add nn 1))
43 (setf np2 (f2cl-lib:int-add nn 2))
44 (setf n2p3 (f2cl-lib:int-add (f2cl-lib:int-mul 2 nn) 3))
45 (setf n3p4 (f2cl-lib:int-add (f2cl-lib:int-mul 3 nn) 4))
46 (setf n4p5 (f2cl-lib:int-add (f2cl-lib:int-mul 4 nn) 5))
47 (setf n5p6 (f2cl-lib:int-add (f2cl-lib:int-mul 5 nn) 6))
48 (setf k (idamax np1 start 1))
49 (setf startk
50 (f2cl-lib:fref start-%data%
51 (k)
52 ((1 (f2cl-lib:int-add nn 1)))
53 start-%offset%))
54 (dcopy lenaa aa 1
55 (f2cl-lib:array-slice work-%data%
56 double-float
57 (n5p6)
58 ((1
59 (f2cl-lib:int-add
60 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
61 lenaa
62 nn
63 1)))
64 work-%offset%)
65 1)
66 (setf (f2cl-lib:fref maxa-%data%
67 ((f2cl-lib:int-add nn 1))
68 ((1 (f2cl-lib:int-add nn 2)))
69 maxa-%offset%)
70 (f2cl-lib:int-add lenaa 1))
71 (setf (f2cl-lib:fref maxa-%data%
72 ((f2cl-lib:int-add nn 2))
73 ((1 (f2cl-lib:int-add nn 2)))
74 maxa-%offset%)
75 (f2cl-lib:int-sub (f2cl-lib:int-add lenaa nn 3) k))
76 (mfacds nn
77 (f2cl-lib:array-slice work-%data%
78 double-float
79 (n5p6)
80 ((1
81 (f2cl-lib:int-add
82 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
83 lenaa
84 nn
85 1)))
86 work-%offset%)
87 lenaa maxa)
88 (dcopy nn pp 1 work 1)
89 (if (< k np1)
90 (setf (f2cl-lib:fref work-%data%
91 (k)
92 ((1
93 (f2cl-lib:int-add
94 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
95 lenaa
96 nn
97 1)))
98 work-%offset%)
99 (+
100 (f2cl-lib:fref work-%data%
101 (k)
102 ((1
103 (f2cl-lib:int-add
104 (f2cl-lib:int-mul 5
105 (f2cl-lib:int-add nn 1))
106 lenaa
107 nn
108 1)))
109 work-%offset%)
110 1.0)))
111 (setf unrm (dnrm2 nn work 1))
112 (setf imax (f2cl-lib:int-mul 10 np1))
113 (setf stillu f2cl-lib:%true%)
114 (setf stillb f2cl-lib:%true%)
115 (setf ztol (* 100.0f0 (f2cl-lib:d1mach 4)))
116 (setf rbtol
117 (* ztol
118 (f2cl-lib:fsqrt
119 (+ (expt startk 2) (expt (dnrm2 nn rho 1) 2)))))
120 (setf rutol (* ztol unrm))
121 (multds
122 (f2cl-lib:array-slice work-%data%
123 double-float
124 (n3p4)
125 ((1
126 (f2cl-lib:int-add
127 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
128 lenaa
129 nn
130 1)))
131 work-%offset%)
132 aa
133 (f2cl-lib:array-slice work-%data%
134 double-float
135 (np2)
136 ((1
137 (f2cl-lib:int-add
138 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
139 lenaa
140 nn
141 1)))
142 work-%offset%)
143 maxa nn lenaa)
144 (setf (f2cl-lib:fref work-%data%
145 ((f2cl-lib:int-add n3p4 nn))
146 ((1
147 (f2cl-lib:int-add
148 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
149 lenaa
150 nn
151 1)))
152 work-%offset%)
153 (f2cl-lib:fref work-%data%
154 ((f2cl-lib:int-sub (f2cl-lib:int-add np2 k) 1))
155 ((1
156 (f2cl-lib:int-add
157 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
158 lenaa
159 nn
160 1)))
161 work-%offset%))
162 (setf ind (f2cl-lib:int-sub (f2cl-lib:int-add n3p4 k) 1))
163 (if (< k np1)
164 (setf (f2cl-lib:fref work-%data%
165 (ind)
166 ((1
167 (f2cl-lib:int-add
168 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
169 lenaa
170 nn
171 1)))
172 work-%offset%)
173 (+
174 (f2cl-lib:fref work-%data%
175 (ind)
176 ((1
177 (f2cl-lib:int-add
178 (f2cl-lib:int-mul 5
179 (f2cl-lib:int-add nn 1))
180 lenaa
181 nn
182 1)))
183 work-%offset%)
184 (f2cl-lib:fref work-%data%
185 ((f2cl-lib:int-add np2 nn))
186 ((1
187 (f2cl-lib:int-add
188 (f2cl-lib:int-mul 5
189 (f2cl-lib:int-add nn 1))
190 lenaa
191 nn
192 1)))
193 work-%offset%))))
194 (dscal np1 -1.0
195 (f2cl-lib:array-slice work-%data%
196 double-float
197 (n3p4)
198 ((1
199 (f2cl-lib:int-add
200 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
201 lenaa
202 nn
203 1)))
204 work-%offset%)
205 1)
206 (daxpy nn -1.0 pp 1
207 (f2cl-lib:array-slice work-%data%
208 double-float
209 (n3p4)
210 ((1
211 (f2cl-lib:int-add
212 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
213 lenaa
214 nn
215 1)))
216 work-%offset%)
217 1)
218 (if (< k np1)
219 (setf (f2cl-lib:fref work-%data%
220 (ind)
221 ((1
222 (f2cl-lib:int-add
223 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
224 lenaa
225 nn
226 1)))
227 work-%offset%)
228 (-
229 (f2cl-lib:fref work-%data%
230 (ind)
231 ((1
232 (f2cl-lib:int-add
233 (f2cl-lib:int-mul 5
234 (f2cl-lib:int-add nn 1))
235 lenaa
236 nn
237 1)))
238 work-%offset%)
239 1.0)))
240 (qimuds
241 (f2cl-lib:array-slice work-%data%
242 double-float
243 (n5p6)
244 ((1
245 (f2cl-lib:int-add
246 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
247 lenaa
248 nn
249 1)))
250 work-%offset%)
251 (f2cl-lib:array-slice work-%data%
252 double-float
253 (n3p4)
254 ((1
255 (f2cl-lib:int-add
256 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
257 lenaa
258 nn
259 1)))
260 work-%offset%)
261 maxa nn lenaa)
262 (dcopy np1
263 (f2cl-lib:array-slice work-%data%
264 double-float
265 (n3p4)
266 ((1
267 (f2cl-lib:int-add
268 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
269 lenaa
270 nn
271 1)))
272 work-%offset%)
273 1 work 1)
274 (qimuds
275 (f2cl-lib:array-slice work-%data%
276 double-float
277 (n5p6)
278 ((1
279 (f2cl-lib:int-add
280 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
281 lenaa
282 nn
283 1)))
284 work-%offset%)
285 work maxa nn lenaa)
286 (multds
287 (f2cl-lib:array-slice work-%data%
288 double-float
289 (n4p5)
290 ((1
291 (f2cl-lib:int-add
292 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
293 lenaa
294 nn
295 1)))
296 work-%offset%)
297 aa work maxa nn lenaa)
298 (setf (f2cl-lib:fref work-%data%
299 ((f2cl-lib:int-add n4p5 nn))
300 ((1
301 (f2cl-lib:int-add
302 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
303 lenaa
304 nn
305 1)))
306 work-%offset%)
307 (f2cl-lib:fref work-%data%
308 (k)
309 ((1
310 (f2cl-lib:int-add
311 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
312 lenaa
313 nn
314 1)))
315 work-%offset%))
316 (if (< k np1)
317 (setf (f2cl-lib:fref work-%data%
318 ((f2cl-lib:int-sub (f2cl-lib:int-add n4p5 k) 1))
319 ((1
320 (f2cl-lib:int-add
321 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
322 lenaa
323 nn
324 1)))
325 work-%offset%)
326 (+
327 (f2cl-lib:fref work-%data%
328 ((f2cl-lib:int-sub (f2cl-lib:int-add n4p5 k)
329 1))
330 ((1
331 (f2cl-lib:int-add
332 (f2cl-lib:int-mul 5
333 (f2cl-lib:int-add nn 1))
334 lenaa
335 nn
336 1)))
337 work-%offset%)
338 (f2cl-lib:fref work-%data%
339 (np1)
340 ((1
341 (f2cl-lib:int-add
342 (f2cl-lib:int-mul 5
343 (f2cl-lib:int-add nn 1))
344 lenaa
345 nn
346 1)))
347 work-%offset%))))
348 (setf runprd
349 (ddot np1
350 (f2cl-lib:array-slice work-%data%
351 double-float
352 (n3p4)
353 ((1
354 (f2cl-lib:int-add
355 (f2cl-lib:int-mul 5
356 (f2cl-lib:int-add nn
357 1))
358 lenaa
359 nn
360 1)))
361 work-%offset%)
362 1
363 (f2cl-lib:array-slice work-%data%
364 double-float
365 (n3p4)
366 ((1
367 (f2cl-lib:int-add
368 (f2cl-lib:int-mul 5
369 (f2cl-lib:int-add nn
370 1))
371 lenaa
372 nn
373 1)))
374 work-%offset%)
375 1))
376 (setf punprd
377 (ddot np1
378 (f2cl-lib:array-slice work-%data%
379 double-float
380 (n4p5)
381 ((1
382 (f2cl-lib:int-add
383 (f2cl-lib:int-mul 5
384 (f2cl-lib:int-add nn
385 1))
386 lenaa
387 nn
388 1)))
389 work-%offset%)
390 1
391 (f2cl-lib:array-slice work-%data%
392 double-float
393 (n4p5)
394 ((1
395 (f2cl-lib:int-add
396 (f2cl-lib:int-mul 5
397 (f2cl-lib:int-add nn
398 1))
399 lenaa
400 nn
401 1)))
402 work-%offset%)
403 1))
404 (setf j 1)
405 label100
406 (if (not (and stillu (<= j imax))) (go label200))
407 (cond
408 ((> (f2cl-lib:fsqrt runprd) rutol)
409 (cond
410 ((= punprd 0.0f0)
411 (multds
412 (f2cl-lib:array-slice work-%data%
413 double-float
414 (n3p4)
415 ((1
416 (f2cl-lib:int-add
417 (f2cl-lib:int-mul 5
418 (f2cl-lib:int-add nn
419 1))
420 lenaa
421 nn
422 1)))
423 work-%offset%)
424 aa
425 (f2cl-lib:array-slice work-%data%
426 double-float
427 (np2)
428 ((1
429 (f2cl-lib:int-add
430 (f2cl-lib:int-mul 5
431 (f2cl-lib:int-add nn
432 1))
433 lenaa
434 nn
435 1)))
436 work-%offset%)
437 maxa nn lenaa)
438 (setf (f2cl-lib:fref work-%data%
439 ((f2cl-lib:int-add n3p4 nn))
440 ((1
441 (f2cl-lib:int-add
442 (f2cl-lib:int-mul 5
443 (f2cl-lib:int-add nn 1))
444 lenaa
445 nn
446 1)))
447 work-%offset%)
448 (f2cl-lib:fref work-%data%
449 ((f2cl-lib:int-sub (f2cl-lib:int-add np2 k)
450 1))
451 ((1
452 (f2cl-lib:int-add
453 (f2cl-lib:int-mul 5
454 (f2cl-lib:int-add nn
455 1))
456 lenaa
457 nn
458 1)))
459 work-%offset%))
460 (setf ind (f2cl-lib:int-sub (f2cl-lib:int-add n3p4 k) 1))
461 (if (< k np1)
462 (setf (f2cl-lib:fref work-%data%
463 (ind)
464 ((1
465 (f2cl-lib:int-add
466 (f2cl-lib:int-mul 5
467 (f2cl-lib:int-add nn
468 1))
469 lenaa
470 nn
471 1)))
472 work-%offset%)
473 (+
474 (f2cl-lib:fref work-%data%
475 (ind)
476 ((1
477 (f2cl-lib:int-add
478 (f2cl-lib:int-mul 5
479 (f2cl-lib:int-add
480 nn
481 1))
482 lenaa
483 nn
484 1)))
485 work-%offset%)
486 (f2cl-lib:fref work-%data%
487 ((f2cl-lib:int-add np2 nn))
488 ((1
489 (f2cl-lib:int-add
490 (f2cl-lib:int-mul 5
491 (f2cl-lib:int-add
492 nn
493 1))
494 lenaa
495 nn
496 1)))
497 work-%offset%))))
498 (dscal np1 -1.0
499 (f2cl-lib:array-slice work-%data%
500 double-float
501 (n3p4)
502 ((1
503 (f2cl-lib:int-add
504 (f2cl-lib:int-mul 5
505 (f2cl-lib:int-add nn
506 1))
507 lenaa
508 nn
509 1)))
510 work-%offset%)
511 1)
512 (daxpy nn -1.0 pp 1
513 (f2cl-lib:array-slice work-%data%
514 double-float
515 (n3p4)
516 ((1
517 (f2cl-lib:int-add
518 (f2cl-lib:int-mul 5
519 (f2cl-lib:int-add nn
520 1))
521 lenaa
522 nn
523 1)))
524 work-%offset%)
525 1)
526 (if (< k np1)
527 (setf (f2cl-lib:fref work-%data%
528 ((f2cl-lib:int-sub
529 (f2cl-lib:int-add n3p4 k)
530 1))
531 ((1
532 (f2cl-lib:int-add
533 (f2cl-lib:int-mul 5
534 (f2cl-lib:int-add nn
535 1))
536 lenaa
537 nn
538 1)))
539 work-%offset%)
540 (-
541 (f2cl-lib:fref work-%data%
542 ((f2cl-lib:int-sub
543 (f2cl-lib:int-add n3p4 k)
544 1))
545 ((1
546 (f2cl-lib:int-add
547 (f2cl-lib:int-mul 5
548 (f2cl-lib:int-add
549 nn
550 1))
551 lenaa
552 nn
553 1)))
554 work-%offset%)
555 1.0)))
556 (qimuds
557 (f2cl-lib:array-slice work-%data%
558 double-float
559 (n5p6)
560 ((1
561 (f2cl-lib:int-add
562 (f2cl-lib:int-mul 5
563 (f2cl-lib:int-add nn
564 1))
565 lenaa
566 nn
567 1)))
568 work-%offset%)
569 (f2cl-lib:array-slice work-%data%
570 double-float
571 (n3p4)
572 ((1
573 (f2cl-lib:int-add
574 (f2cl-lib:int-mul 5
575 (f2cl-lib:int-add nn
576 1))
577 lenaa
578 nn
579 1)))
580 work-%offset%)
581 maxa nn lenaa)
582 (dcopy np1
583 (f2cl-lib:array-slice work-%data%
584 double-float
585 (n3p4)
586 ((1
587 (f2cl-lib:int-add
588 (f2cl-lib:int-mul 5
589 (f2cl-lib:int-add nn
590 1))
591 lenaa
592 nn
593 1)))
594 work-%offset%)
595 1 work 1)
596 (qimuds
597 (f2cl-lib:array-slice work-%data%
598 double-float
599 (n5p6)
600 ((1
601 (f2cl-lib:int-add
602 (f2cl-lib:int-mul 5
603 (f2cl-lib:int-add nn
604 1))
605 lenaa
606 nn
607 1)))
608 work-%offset%)
609 work maxa nn lenaa)
610 (multds
611 (f2cl-lib:array-slice work-%data%
612 double-float
613 (n4p5)
614 ((1
615 (f2cl-lib:int-add
616 (f2cl-lib:int-mul 5
617 (f2cl-lib:int-add nn
618 1))
619 lenaa
620 nn
621 1)))
622 work-%offset%)
623 aa work maxa nn lenaa)
624 (setf (f2cl-lib:fref work-%data%
625 ((f2cl-lib:int-add n4p5 nn))
626 ((1
627 (f2cl-lib:int-add
628 (f2cl-lib:int-mul 5
629 (f2cl-lib:int-add nn 1))
630 lenaa
631 nn
632 1)))
633 work-%offset%)
634 (f2cl-lib:fref work-%data%
635 (k)
636 ((1
637 (f2cl-lib:int-add
638 (f2cl-lib:int-mul 5
639 (f2cl-lib:int-add nn
640 1))
641 lenaa
642 nn
643 1)))
644 work-%offset%))
645 (setf ind (f2cl-lib:int-sub (f2cl-lib:int-add n4p5 k) 1))
646 (if (< k np1)
647 (setf (f2cl-lib:fref work-%data%
648 (ind)
649 ((1
650 (f2cl-lib:int-add
651 (f2cl-lib:int-mul 5
652 (f2cl-lib:int-add nn
653 1))
654 lenaa
655 nn
656 1)))
657 work-%offset%)
658 (+
659 (f2cl-lib:fref work-%data%
660 (ind)
661 ((1
662 (f2cl-lib:int-add
663 (f2cl-lib:int-mul 5
664 (f2cl-lib:int-add
665 nn
666 1))
667 lenaa
668 nn
669 1)))
670 work-%offset%)
671 (f2cl-lib:fref work-%data%
672 (np1)
673 ((1
674 (f2cl-lib:int-add
675 (f2cl-lib:int-mul 5
676 (f2cl-lib:int-add
677 nn
678 1))
679 lenaa
680 nn
681 1)))
682 work-%offset%))))
683 (setf runprd
684 (ddot np1
685 (f2cl-lib:array-slice work-%data%
686 double-float
687 (n3p4)
688 ((1
689 (f2cl-lib:int-add
690 (f2cl-lib:int-mul 5
691 (f2cl-lib:int-add
692 nn
693 1))
694 lenaa
695 nn
696 1)))
697 work-%offset%)
698 1
699 (f2cl-lib:array-slice work-%data%
700 double-float
701 (n3p4)
702 ((1
703 (f2cl-lib:int-add
704 (f2cl-lib:int-mul 5
705 (f2cl-lib:int-add
706 nn
707 1))
708 lenaa
709 nn
710 1)))
711 work-%offset%)
712 1))
713 (setf punprd
714 (ddot np1
715 (f2cl-lib:array-slice work-%data%
716 double-float
717 (n4p5)
718 ((1
719 (f2cl-lib:int-add
720 (f2cl-lib:int-mul 5
721 (f2cl-lib:int-add
722 nn
723 1))
724 lenaa
725 nn
726 1)))
727 work-%offset%)
728 1
729 (f2cl-lib:array-slice work-%data%
730 double-float
731 (n4p5)
732 ((1
733 (f2cl-lib:int-add
734 (f2cl-lib:int-mul 5
735 (f2cl-lib:int-add
736 nn
737 1))
738 lenaa
739 nn
740 1)))
741 work-%offset%)
742 1))
743 (cond
744 ((<= (f2cl-lib:fsqrt runprd) rutol)
745 (setf stillu f2cl-lib:%false%)))))
746 (cond
747 (stillu
748 (setf au (/ runprd punprd))
749 (dcopy np1
750 (f2cl-lib:array-slice work-%data%
751 double-float
752 (np2)
753 ((1
754 (f2cl-lib:int-add
755 (f2cl-lib:int-mul 5
756 (f2cl-lib:int-add nn
757 1))
758 lenaa
759 nn
760 1)))
761 work-%offset%)
762 1 work 1)
763 (daxpy np1 au
764 (f2cl-lib:array-slice work-%data%
765 double-float
766 (n4p5)
767 ((1
768 (f2cl-lib:int-add
769 (f2cl-lib:int-mul 5
770 (f2cl-lib:int-add nn
771 1))
772 lenaa
773 nn
774 1)))
775 work-%offset%)
776 1
777 (f2cl-lib:array-slice work-%data%
778 double-float
779 (np2)
780 ((1
781 (f2cl-lib:int-add
782 (f2cl-lib:int-mul 5
783 (f2cl-lib:int-add nn
784 1))
785 lenaa
786 nn
787 1)))
788 work-%offset%)
789 1)
790 (daxpy np1 -1.0
791 (f2cl-lib:array-slice work-%data%
792 double-float
793 (np2)
794 ((1
795 (f2cl-lib:int-add
796 (f2cl-lib:int-mul 5
797 (f2cl-lib:int-add nn
798 1))
799 lenaa
800 nn
801 1)))
802 work-%offset%)
803 1 work 1)
804 (setf zlen
805 (dnrm2 np1
806 (f2cl-lib:array-slice work-%data%
807 double-float
808 (np2)
809 ((1
810 (f2cl-lib:int-add
811 (f2cl-lib:int-mul 5
812 (f2cl-lib:int-add
813 nn
814 1))
815 lenaa
816 nn
817 1)))
818 work-%offset%)
819 1))
820 (setf dznrm (dnrm2 np1 work 1))
821 (if (< (/ dznrm zlen) ztol) (setf stillu f2cl-lib:%false%)))))
822 (t
823 (setf stillu f2cl-lib:%false%)))
824 (cond
825 (stillu
826 (multds work aa
827 (f2cl-lib:array-slice work-%data%
828 double-float
829 (n4p5)
830 ((1
831 (f2cl-lib:int-add
832 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
833 lenaa
834 nn
835 1)))
836 work-%offset%)
837 maxa nn lenaa)
838 (setf (f2cl-lib:fref work-%data%
839 (np1)
840 ((1
841 (f2cl-lib:int-add
842 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
843 lenaa
844 nn
845 1)))
846 work-%offset%)
847 (f2cl-lib:fref work-%data%
848 ((f2cl-lib:int-sub (f2cl-lib:int-add n4p5 k)
849 1))
850 ((1
851 (f2cl-lib:int-add
852 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
853 lenaa
854 nn
855 1)))
856 work-%offset%))
857 (if (< k np1)
858 (setf (f2cl-lib:fref work-%data%
859 (k)
860 ((1
861 (f2cl-lib:int-add
862 (f2cl-lib:int-mul 5
863 (f2cl-lib:int-add nn 1))
864 lenaa
865 nn
866 1)))
867 work-%offset%)
868 (+
869 (f2cl-lib:fref work-%data%
870 (k)
871 ((1
872 (f2cl-lib:int-add
873 (f2cl-lib:int-mul 5
874 (f2cl-lib:int-add nn
875 1))
876 lenaa
877 nn
878 1)))
879 work-%offset%)
880 (f2cl-lib:fref work-%data%
881 ((f2cl-lib:int-add n4p5 nn))
882 ((1
883 (f2cl-lib:int-add
884 (f2cl-lib:int-mul 5
885 (f2cl-lib:int-add nn
886 1))
887 lenaa
888 nn
889 1)))
890 work-%offset%))))
891 (qimuds
892 (f2cl-lib:array-slice work-%data%
893 double-float
894 (n5p6)
895 ((1
896 (f2cl-lib:int-add
897 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
898 lenaa
899 nn
900 1)))
901 work-%offset%)
902 work maxa nn lenaa)
903 (daxpy np1 (- au) work 1
904 (f2cl-lib:array-slice work-%data%
905 double-float
906 (n3p4)
907 ((1
908 (f2cl-lib:int-add
909 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
910 lenaa
911 nn
912 1)))
913 work-%offset%)
914 1)
915 (setf rnprd
916 (ddot np1
917 (f2cl-lib:array-slice work-%data%
918 double-float
919 (n3p4)
920 ((1
921 (f2cl-lib:int-add
922 (f2cl-lib:int-mul 5
923 (f2cl-lib:int-add
924 nn
925 1))
926 lenaa
927 nn
928 1)))
929 work-%offset%)
930 1
931 (f2cl-lib:array-slice work-%data%
932 double-float
933 (n3p4)
934 ((1
935 (f2cl-lib:int-add
936 (f2cl-lib:int-mul 5
937 (f2cl-lib:int-add
938 nn
939 1))
940 lenaa
941 nn
942 1)))
943 work-%offset%)
944 1))
945 (setf bu (/ rnprd runprd))
946 (setf runprd rnprd)
947 (dcopy np1
948 (f2cl-lib:array-slice work-%data%
949 double-float
950 (n3p4)
951 ((1
952 (f2cl-lib:int-add
953 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
954 lenaa
955 nn
956 1)))
957 work-%offset%)
958 1 work 1)
959 (qimuds
960 (f2cl-lib:array-slice work-%data%
961 double-float
962 (n5p6)
963 ((1
964 (f2cl-lib:int-add
965 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
966 lenaa
967 nn
968 1)))
969 work-%offset%)
970 work maxa nn lenaa)
971 (multds start aa work maxa nn lenaa)
972 (setf (f2cl-lib:fref start-%data%
973 (np1)
974 ((1 (f2cl-lib:int-add nn 1)))
975 start-%offset%)
976 (f2cl-lib:fref work-%data%
977 (k)
978 ((1
979 (f2cl-lib:int-add
980 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
981 lenaa
982 nn
983 1)))
984 work-%offset%))
985 (if (< k np1)
986 (setf (f2cl-lib:fref start-%data%
987 (k)
988 ((1 (f2cl-lib:int-add nn 1)))
989 start-%offset%)
990 (+
991 (f2cl-lib:fref start-%data%
992 (k)
993 ((1 (f2cl-lib:int-add nn 1)))
994 start-%offset%)
995 (f2cl-lib:fref work-%data%
996 (np1)
997 ((1
998 (f2cl-lib:int-add
999 (f2cl-lib:int-mul 5
1000 (f2cl-lib:int-add nn
1001 1))
1002 lenaa
1003 nn
1004 1)))
1005 work-%offset%))))
1006 (daxpy np1 bu
1007 (f2cl-lib:array-slice work-%data%
1008 double-float
1009 (n4p5)
1010 ((1
1011 (f2cl-lib:int-add
1012 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1013 lenaa
1014 nn
1015 1)))
1016 work-%offset%)
1017 1 start 1)
1018 (dcopy np1 start 1
1019 (f2cl-lib:array-slice work-%data%
1020 double-float
1021 (n4p5)
1022 ((1
1023 (f2cl-lib:int-add
1024 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1025 lenaa
1026 nn
1027 1)))
1028 work-%offset%)
1029 1)
1030 (setf punprd
1031 (ddot np1
1032 (f2cl-lib:array-slice work-%data%
1033 double-float
1034 (n4p5)
1035 ((1
1036 (f2cl-lib:int-add
1037 (f2cl-lib:int-mul 5
1038 (f2cl-lib:int-add
1039 nn
1040 1))
1041 lenaa
1042 nn
1043 1)))
1044 work-%offset%)
1045 1
1046 (f2cl-lib:array-slice work-%data%
1047 double-float
1048 (n4p5)
1049 ((1
1050 (f2cl-lib:int-add
1051 (f2cl-lib:int-mul 5
1052 (f2cl-lib:int-add
1053 nn
1054 1))
1055 lenaa
1056 nn
1057 1)))
1058 work-%offset%)
1059 1))))
1060 (setf j (f2cl-lib:int-add j 1))
1061 (go label100)
1062 label200
1063 (cond
1064 ((> j imax)
1065 (setf iflag 4)
1066 (go end_label)))
1067 (multds
1068 (f2cl-lib:array-slice work-%data%
1069 double-float
1070 (n3p4)
1071 ((1
1072 (f2cl-lib:int-add
1073 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1074 lenaa
1075 nn
1076 1)))
1077 work-%offset%)
1078 aa
1079 (f2cl-lib:array-slice work-%data%
1080 double-float
1081 (n2p3)
1082 ((1
1083 (f2cl-lib:int-add
1084 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1085 lenaa
1086 nn
1087 1)))
1088 work-%offset%)
1089 maxa nn lenaa)
1090 (setf (f2cl-lib:fref work-%data%
1091 ((f2cl-lib:int-add n3p4 nn))
1092 ((1
1093 (f2cl-lib:int-add
1094 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1095 lenaa
1096 nn
1097 1)))
1098 work-%offset%)
1099 (f2cl-lib:fref work-%data%
1100 ((f2cl-lib:int-sub (f2cl-lib:int-add n2p3 k) 1))
1101 ((1
1102 (f2cl-lib:int-add
1103 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1104 lenaa
1105 nn
1106 1)))
1107 work-%offset%))
1108 (setf ind (f2cl-lib:int-sub (f2cl-lib:int-add n3p4 k) 1))
1109 (if (< k np1)
1110 (setf (f2cl-lib:fref work-%data%
1111 (ind)
1112 ((1
1113 (f2cl-lib:int-add
1114 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1115 lenaa
1116 nn
1117 1)))
1118 work-%offset%)
1119 (+
1120 (f2cl-lib:fref work-%data%
1121 (ind)
1122 ((1
1123 (f2cl-lib:int-add
1124 (f2cl-lib:int-mul 5
1125 (f2cl-lib:int-add nn 1))
1126 lenaa
1127 nn
1128 1)))
1129 work-%offset%)
1130 (f2cl-lib:fref work-%data%
1131 ((f2cl-lib:int-add n2p3 nn))
1132 ((1
1133 (f2cl-lib:int-add
1134 (f2cl-lib:int-mul 5
1135 (f2cl-lib:int-add nn 1))
1136 lenaa
1137 nn
1138 1)))
1139 work-%offset%))))
1140 (dscal np1 -1.0
1141 (f2cl-lib:array-slice work-%data%
1142 double-float
1143 (n3p4)
1144 ((1
1145 (f2cl-lib:int-add
1146 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1147 lenaa
1148 nn
1149 1)))
1150 work-%offset%)
1151 1)
1152 (daxpy nn -1.0 rho 1
1153 (f2cl-lib:array-slice work-%data%
1154 double-float
1155 (n3p4)
1156 ((1
1157 (f2cl-lib:int-add
1158 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1159 lenaa
1160 nn
1161 1)))
1162 work-%offset%)
1163 1)
1164 (setf (f2cl-lib:fref work-%data%
1165 ((f2cl-lib:int-add n3p4 nn))
1166 ((1
1167 (f2cl-lib:int-add
1168 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1169 lenaa
1170 nn
1171 1)))
1172 work-%offset%)
1173 (+ startk
1174 (f2cl-lib:fref work-%data%
1175 ((f2cl-lib:int-add n3p4 nn))
1176 ((1
1177 (f2cl-lib:int-add
1178 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1179 lenaa
1180 nn
1181 1)))
1182 work-%offset%)))
1183 (qimuds
1184 (f2cl-lib:array-slice work-%data%
1185 double-float
1186 (n5p6)
1187 ((1
1188 (f2cl-lib:int-add
1189 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1190 lenaa
1191 nn
1192 1)))
1193 work-%offset%)
1194 (f2cl-lib:array-slice work-%data%
1195 double-float
1196 (n3p4)
1197 ((1
1198 (f2cl-lib:int-add
1199 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1200 lenaa
1201 nn
1202 1)))
1203 work-%offset%)
1204 maxa nn lenaa)
1205 (dcopy np1
1206 (f2cl-lib:array-slice work-%data%
1207 double-float
1208 (n3p4)
1209 ((1
1210 (f2cl-lib:int-add
1211 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1212 lenaa
1213 nn
1214 1)))
1215 work-%offset%)
1216 1 work 1)
1217 (qimuds
1218 (f2cl-lib:array-slice work-%data%
1219 double-float
1220 (n5p6)
1221 ((1
1222 (f2cl-lib:int-add
1223 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1224 lenaa
1225 nn
1226 1)))
1227 work-%offset%)
1228 work maxa nn lenaa)
1229 (multds
1230 (f2cl-lib:array-slice work-%data%
1231 double-float
1232 (n4p5)
1233 ((1
1234 (f2cl-lib:int-add
1235 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1236 lenaa
1237 nn
1238 1)))
1239 work-%offset%)
1240 aa work maxa nn lenaa)
1241 (setf (f2cl-lib:fref work-%data%
1242 ((f2cl-lib:int-add n4p5 nn))
1243 ((1
1244 (f2cl-lib:int-add
1245 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1246 lenaa
1247 nn
1248 1)))
1249 work-%offset%)
1250 (f2cl-lib:fref work-%data%
1251 (k)
1252 ((1
1253 (f2cl-lib:int-add
1254 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1255 lenaa
1256 nn
1257 1)))
1258 work-%offset%))
1259 (if (< k np1)
1260 (setf (f2cl-lib:fref work-%data%
1261 ((f2cl-lib:int-sub (f2cl-lib:int-add n4p5 k) 1))
1262 ((1
1263 (f2cl-lib:int-add
1264 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1265 lenaa
1266 nn
1267 1)))
1268 work-%offset%)
1269 (+
1270 (f2cl-lib:fref work-%data%
1271 ((f2cl-lib:int-sub (f2cl-lib:int-add n4p5 k)
1272 1))
1273 ((1
1274 (f2cl-lib:int-add
1275 (f2cl-lib:int-mul 5
1276 (f2cl-lib:int-add nn 1))
1277 lenaa
1278 nn
1279 1)))
1280 work-%offset%)
1281 (f2cl-lib:fref work-%data%
1282 (np1)
1283 ((1
1284 (f2cl-lib:int-add
1285 (f2cl-lib:int-mul 5
1286 (f2cl-lib:int-add nn 1))
1287 lenaa
1288 nn
1289 1)))
1290 work-%offset%))))
1291 (setf rbnprd
1292 (ddot np1
1293 (f2cl-lib:array-slice work-%data%
1294 double-float
1295 (n3p4)
1296 ((1
1297 (f2cl-lib:int-add
1298 (f2cl-lib:int-mul 5
1299 (f2cl-lib:int-add nn
1300 1))
1301 lenaa
1302 nn
1303 1)))
1304 work-%offset%)
1305 1
1306 (f2cl-lib:array-slice work-%data%
1307 double-float
1308 (n3p4)
1309 ((1
1310 (f2cl-lib:int-add
1311 (f2cl-lib:int-mul 5
1312 (f2cl-lib:int-add nn
1313 1))
1314 lenaa
1315 nn
1316 1)))
1317 work-%offset%)
1318 1))
1319 (setf pbnprd
1320 (ddot np1
1321 (f2cl-lib:array-slice work-%data%
1322 double-float
1323 (n4p5)
1324 ((1
1325 (f2cl-lib:int-add
1326 (f2cl-lib:int-mul 5
1327 (f2cl-lib:int-add nn
1328 1))
1329 lenaa
1330 nn
1331 1)))
1332 work-%offset%)
1333 1
1334 (f2cl-lib:array-slice work-%data%
1335 double-float
1336 (n4p5)
1337 ((1
1338 (f2cl-lib:int-add
1339 (f2cl-lib:int-mul 5
1340 (f2cl-lib:int-add nn
1341 1))
1342 lenaa
1343 nn
1344 1)))
1345 work-%offset%)
1346 1))
1347 (setf j 1)
1348 label300
1349 (if (not (and stillb (<= j imax))) (go label400))
1350 (cond
1351 ((> (f2cl-lib:fsqrt rbnprd) rbtol)
1352 (cond
1353 ((= pbnprd 0.0f0)
1354 (multds
1355 (f2cl-lib:array-slice work-%data%
1356 double-float
1357 (n3p4)
1358 ((1
1359 (f2cl-lib:int-add
1360 (f2cl-lib:int-mul 5
1361 (f2cl-lib:int-add nn
1362 1))
1363 lenaa
1364 nn
1365 1)))
1366 work-%offset%)
1367 aa
1368 (f2cl-lib:array-slice work-%data%
1369 double-float
1370 (n2p3)
1371 ((1
1372 (f2cl-lib:int-add
1373 (f2cl-lib:int-mul 5
1374 (f2cl-lib:int-add nn
1375 1))
1376 lenaa
1377 nn
1378 1)))
1379 work-%offset%)
1380 maxa nn lenaa)
1381 (setf (f2cl-lib:fref work-%data%
1382 ((f2cl-lib:int-add n3p4 nn))
1383 ((1
1384 (f2cl-lib:int-add
1385 (f2cl-lib:int-mul 5
1386 (f2cl-lib:int-add nn 1))
1387 lenaa
1388 nn
1389 1)))
1390 work-%offset%)
1391 (f2cl-lib:fref work-%data%
1392 ((f2cl-lib:int-sub (f2cl-lib:int-add n2p3 k)
1393 1))
1394 ((1
1395 (f2cl-lib:int-add
1396 (f2cl-lib:int-mul 5
1397 (f2cl-lib:int-add nn
1398 1))
1399 lenaa
1400 nn
1401 1)))
1402 work-%offset%))
1403 (setf ind (f2cl-lib:int-sub (f2cl-lib:int-add n3p4 k) 1))
1404 (if (< k np1)
1405 (setf (f2cl-lib:fref work-%data%
1406 (ind)
1407 ((1
1408 (f2cl-lib:int-add
1409 (f2cl-lib:int-mul 5
1410 (f2cl-lib:int-add nn
1411 1))
1412 lenaa
1413 nn
1414 1)))
1415 work-%offset%)
1416 (+
1417 (f2cl-lib:fref work-%data%
1418 (ind)
1419 ((1
1420 (f2cl-lib:int-add
1421 (f2cl-lib:int-mul 5
1422 (f2cl-lib:int-add
1423 nn
1424 1))
1425 lenaa
1426 nn
1427 1)))
1428 work-%offset%)
1429 (f2cl-lib:fref work-%data%
1430 ((f2cl-lib:int-add n2p3 nn))
1431 ((1
1432 (f2cl-lib:int-add
1433 (f2cl-lib:int-mul 5
1434 (f2cl-lib:int-add
1435 nn
1436 1))
1437 lenaa
1438 nn
1439 1)))
1440 work-%offset%))))
1441 (dscal np1 -1.0
1442 (f2cl-lib:array-slice work-%data%
1443 double-float
1444 (n3p4)
1445 ((1
1446 (f2cl-lib:int-add
1447 (f2cl-lib:int-mul 5
1448 (f2cl-lib:int-add nn
1449 1))
1450 lenaa
1451 nn
1452 1)))
1453 work-%offset%)
1454 1)
1455 (daxpy nn -1.0 rho 1
1456 (f2cl-lib:array-slice work-%data%
1457 double-float
1458 (n3p4)
1459 ((1
1460 (f2cl-lib:int-add
1461 (f2cl-lib:int-mul 5
1462 (f2cl-lib:int-add nn
1463 1))
1464 lenaa
1465 nn
1466 1)))
1467 work-%offset%)
1468 1)
1469 (setf (f2cl-lib:fref work-%data%
1470 ((f2cl-lib:int-add n3p4 nn))
1471 ((1
1472 (f2cl-lib:int-add
1473 (f2cl-lib:int-mul 5
1474 (f2cl-lib:int-add nn 1))
1475 lenaa
1476 nn
1477 1)))
1478 work-%offset%)
1479 (+ startk
1480 (f2cl-lib:fref work-%data%
1481 ((f2cl-lib:int-add n3p4 nn))
1482 ((1
1483 (f2cl-lib:int-add
1484 (f2cl-lib:int-mul 5
1485 (f2cl-lib:int-add nn
1486 1))
1487 lenaa
1488 nn
1489 1)))
1490 work-%offset%)))
1491 (qimuds
1492 (f2cl-lib:array-slice work-%data%
1493 double-float
1494 (n5p6)
1495 ((1
1496 (f2cl-lib:int-add
1497 (f2cl-lib:int-mul 5
1498 (f2cl-lib:int-add nn
1499 1))
1500 lenaa
1501 nn
1502 1)))
1503 work-%offset%)
1504 (f2cl-lib:array-slice work-%data%
1505 double-float
1506 (n3p4)
1507 ((1
1508 (f2cl-lib:int-add
1509 (f2cl-lib:int-mul 5
1510 (f2cl-lib:int-add nn
1511 1))
1512 lenaa
1513 nn
1514 1)))
1515 work-%offset%)
1516 maxa nn lenaa)
1517 (dcopy np1
1518 (f2cl-lib:array-slice work-%data%
1519 double-float
1520 (n3p4)
1521 ((1
1522 (f2cl-lib:int-add
1523 (f2cl-lib:int-mul 5
1524 (f2cl-lib:int-add nn
1525 1))
1526 lenaa
1527 nn
1528 1)))
1529 work-%offset%)
1530 1 work 1)
1531 (qimuds
1532 (f2cl-lib:array-slice work-%data%
1533 double-float
1534 (n5p6)
1535 ((1
1536 (f2cl-lib:int-add
1537 (f2cl-lib:int-mul 5
1538 (f2cl-lib:int-add nn
1539 1))
1540 lenaa
1541 nn
1542 1)))
1543 work-%offset%)
1544 work maxa nn lenaa)
1545 (multds
1546 (f2cl-lib:array-slice work-%data%
1547 double-float
1548 (n4p5)
1549 ((1
1550 (f2cl-lib:int-add
1551 (f2cl-lib:int-mul 5
1552 (f2cl-lib:int-add nn
1553 1))
1554 lenaa
1555 nn
1556 1)))
1557 work-%offset%)
1558 aa work maxa nn lenaa)
1559 (setf (f2cl-lib:fref work-%data%
1560 ((f2cl-lib:int-add n4p5 nn))
1561 ((1
1562 (f2cl-lib:int-add
1563 (f2cl-lib:int-mul 5
1564 (f2cl-lib:int-add nn 1))
1565 lenaa
1566 nn
1567 1)))
1568 work-%offset%)
1569 (f2cl-lib:fref work-%data%
1570 (k)
1571 ((1
1572 (f2cl-lib:int-add
1573 (f2cl-lib:int-mul 5
1574 (f2cl-lib:int-add nn
1575 1))
1576 lenaa
1577 nn
1578 1)))
1579 work-%offset%))
1580 (setf ind (f2cl-lib:int-sub (f2cl-lib:int-add n4p5 k) 1))
1581 (if (< k np1)
1582 (setf (f2cl-lib:fref work-%data%
1583 (ind)
1584 ((1
1585 (f2cl-lib:int-add
1586 (f2cl-lib:int-mul 5
1587 (f2cl-lib:int-add nn
1588 1))
1589 lenaa
1590 nn
1591 1)))
1592 work-%offset%)
1593 (+
1594 (f2cl-lib:fref work-%data%
1595 (ind)
1596 ((1
1597 (f2cl-lib:int-add
1598 (f2cl-lib:int-mul 5
1599 (f2cl-lib:int-add
1600 nn
1601 1))
1602 lenaa
1603 nn
1604 1)))
1605 work-%offset%)
1606 (f2cl-lib:fref work-%data%
1607 (np1)
1608 ((1
1609 (f2cl-lib:int-add
1610 (f2cl-lib:int-mul 5
1611 (f2cl-lib:int-add
1612 nn
1613 1))
1614 lenaa
1615 nn
1616 1)))
1617 work-%offset%))))
1618 (setf rbnprd
1619 (ddot np1
1620 (f2cl-lib:array-slice work-%data%
1621 double-float
1622 (n3p4)
1623 ((1
1624 (f2cl-lib:int-add
1625 (f2cl-lib:int-mul 5
1626 (f2cl-lib:int-add
1627 nn
1628 1))
1629 lenaa
1630 nn
1631 1)))
1632 work-%offset%)
1633 1
1634 (f2cl-lib:array-slice work-%data%
1635 double-float
1636 (n3p4)
1637 ((1
1638 (f2cl-lib:int-add
1639 (f2cl-lib:int-mul 5
1640 (f2cl-lib:int-add
1641 nn
1642 1))
1643 lenaa
1644 nn
1645 1)))
1646 work-%offset%)
1647 1))
1648 (setf pbnprd
1649 (ddot np1
1650 (f2cl-lib:array-slice work-%data%
1651 double-float
1652 (n4p5)
1653 ((1
1654 (f2cl-lib:int-add
1655 (f2cl-lib:int-mul 5
1656 (f2cl-lib:int-add
1657 nn
1658 1))
1659 lenaa
1660 nn
1661 1)))
1662 work-%offset%)
1663 1
1664 (f2cl-lib:array-slice work-%data%
1665 double-float
1666 (n4p5)
1667 ((1
1668 (f2cl-lib:int-add
1669 (f2cl-lib:int-mul 5
1670 (f2cl-lib:int-add
1671 nn
1672 1))
1673 lenaa
1674 nn
1675 1)))
1676 work-%offset%)
1677 1))
1678 (cond
1679 ((<= (f2cl-lib:fsqrt rbnprd) rbtol)
1680 (setf stillb f2cl-lib:%false%)))))
1681 (cond
1682 (stillb
1683 (setf ab (/ rbnprd pbnprd))
1684 (dcopy np1
1685 (f2cl-lib:array-slice work-%data%
1686 double-float
1687 (n2p3)
1688 ((1
1689 (f2cl-lib:int-add
1690 (f2cl-lib:int-mul 5
1691 (f2cl-lib:int-add nn
1692 1))
1693 lenaa
1694 nn
1695 1)))
1696 work-%offset%)
1697 1 work 1)
1698 (daxpy np1 ab
1699 (f2cl-lib:array-slice work-%data%
1700 double-float
1701 (n4p5)
1702 ((1
1703 (f2cl-lib:int-add
1704 (f2cl-lib:int-mul 5
1705 (f2cl-lib:int-add nn
1706 1))
1707 lenaa
1708 nn
1709 1)))
1710 work-%offset%)
1711 1
1712 (f2cl-lib:array-slice work-%data%
1713 double-float
1714 (n2p3)
1715 ((1
1716 (f2cl-lib:int-add
1717 (f2cl-lib:int-mul 5
1718 (f2cl-lib:int-add nn
1719 1))
1720 lenaa
1721 nn
1722 1)))
1723 work-%offset%)
1724 1)
1725 (daxpy np1 -1.0
1726 (f2cl-lib:array-slice work-%data%
1727 double-float
1728 (n2p3)
1729 ((1
1730 (f2cl-lib:int-add
1731 (f2cl-lib:int-mul 5
1732 (f2cl-lib:int-add nn
1733 1))
1734 lenaa
1735 nn
1736 1)))
1737 work-%offset%)
1738 1 work 1)
1739 (setf zlen
1740 (dnrm2 np1
1741 (f2cl-lib:array-slice work-%data%
1742 double-float
1743 (n2p3)
1744 ((1
1745 (f2cl-lib:int-add
1746 (f2cl-lib:int-mul 5
1747 (f2cl-lib:int-add
1748 nn
1749 1))
1750 lenaa
1751 nn
1752 1)))
1753 work-%offset%)
1754 1))
1755 (setf dznrm (dnrm2 np1 work 1))
1756 (if (< (/ dznrm zlen) ztol) (setf stillb f2cl-lib:%false%)))))
1757 (t
1758 (setf stillb f2cl-lib:%false%)))
1759 (cond
1760 (stillb
1761 (multds work aa
1762 (f2cl-lib:array-slice work-%data%
1763 double-float
1764 (n4p5)
1765 ((1
1766 (f2cl-lib:int-add
1767 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1768 lenaa
1769 nn
1770 1)))
1771 work-%offset%)
1772 maxa nn lenaa)
1773 (setf (f2cl-lib:fref work-%data%
1774 (np1)
1775 ((1
1776 (f2cl-lib:int-add
1777 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1778 lenaa
1779 nn
1780 1)))
1781 work-%offset%)
1782 (f2cl-lib:fref work-%data%
1783 ((f2cl-lib:int-sub (f2cl-lib:int-add n4p5 k)
1784 1))
1785 ((1
1786 (f2cl-lib:int-add
1787 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1788 lenaa
1789 nn
1790 1)))
1791 work-%offset%))
1792 (if (< k np1)
1793 (setf (f2cl-lib:fref work-%data%
1794 (k)
1795 ((1
1796 (f2cl-lib:int-add
1797 (f2cl-lib:int-mul 5
1798 (f2cl-lib:int-add nn 1))
1799 lenaa
1800 nn
1801 1)))
1802 work-%offset%)
1803 (+
1804 (f2cl-lib:fref work-%data%
1805 (k)
1806 ((1
1807 (f2cl-lib:int-add
1808 (f2cl-lib:int-mul 5
1809 (f2cl-lib:int-add nn
1810 1))
1811 lenaa
1812 nn
1813 1)))
1814 work-%offset%)
1815 (f2cl-lib:fref work-%data%
1816 ((f2cl-lib:int-add n4p5 nn))
1817 ((1
1818 (f2cl-lib:int-add
1819 (f2cl-lib:int-mul 5
1820 (f2cl-lib:int-add nn
1821 1))
1822 lenaa
1823 nn
1824 1)))
1825 work-%offset%))))
1826 (qimuds
1827 (f2cl-lib:array-slice work-%data%
1828 double-float
1829 (n5p6)
1830 ((1
1831 (f2cl-lib:int-add
1832 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1833 lenaa
1834 nn
1835 1)))
1836 work-%offset%)
1837 work maxa nn lenaa)
1838 (daxpy np1 (- ab) work 1
1839 (f2cl-lib:array-slice work-%data%
1840 double-float
1841 (n3p4)
1842 ((1
1843 (f2cl-lib:int-add
1844 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1845 lenaa
1846 nn
1847 1)))
1848 work-%offset%)
1849 1)
1850 (setf rnprd
1851 (ddot np1
1852 (f2cl-lib:array-slice work-%data%
1853 double-float
1854 (n3p4)
1855 ((1
1856 (f2cl-lib:int-add
1857 (f2cl-lib:int-mul 5
1858 (f2cl-lib:int-add
1859 nn
1860 1))
1861 lenaa
1862 nn
1863 1)))
1864 work-%offset%)
1865 1
1866 (f2cl-lib:array-slice work-%data%
1867 double-float
1868 (n3p4)
1869 ((1
1870 (f2cl-lib:int-add
1871 (f2cl-lib:int-mul 5
1872 (f2cl-lib:int-add
1873 nn
1874 1))
1875 lenaa
1876 nn
1877 1)))
1878 work-%offset%)
1879 1))
1880 (setf bb (/ rnprd rbnprd))
1881 (setf rbnprd rnprd)
1882 (dcopy np1
1883 (f2cl-lib:array-slice work-%data%
1884 double-float
1885 (n3p4)
1886 ((1
1887 (f2cl-lib:int-add
1888 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1889 lenaa
1890 nn
1891 1)))
1892 work-%offset%)
1893 1 work 1)
1894 (qimuds
1895 (f2cl-lib:array-slice work-%data%
1896 double-float
1897 (n5p6)
1898 ((1
1899 (f2cl-lib:int-add
1900 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1901 lenaa
1902 nn
1903 1)))
1904 work-%offset%)
1905 work maxa nn lenaa)
1906 (multds start aa work maxa nn lenaa)
1907 (setf (f2cl-lib:fref start-%data%
1908 (np1)
1909 ((1 (f2cl-lib:int-add nn 1)))
1910 start-%offset%)
1911 (f2cl-lib:fref work-%data%
1912 (k)
1913 ((1
1914 (f2cl-lib:int-add
1915 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1916 lenaa
1917 nn
1918 1)))
1919 work-%offset%))
1920 (if (< k np1)
1921 (setf (f2cl-lib:fref start-%data%
1922 (k)
1923 ((1 (f2cl-lib:int-add nn 1)))
1924 start-%offset%)
1925 (+
1926 (f2cl-lib:fref start-%data%
1927 (k)
1928 ((1 (f2cl-lib:int-add nn 1)))
1929 start-%offset%)
1930 (f2cl-lib:fref work-%data%
1931 (np1)
1932 ((1
1933 (f2cl-lib:int-add
1934 (f2cl-lib:int-mul 5
1935 (f2cl-lib:int-add nn
1936 1))
1937 lenaa
1938 nn
1939 1)))
1940 work-%offset%))))
1941 (daxpy np1 bb
1942 (f2cl-lib:array-slice work-%data%
1943 double-float
1944 (n4p5)
1945 ((1
1946 (f2cl-lib:int-add
1947 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1948 lenaa
1949 nn
1950 1)))
1951 work-%offset%)
1952 1 start 1)
1953 (dcopy np1 start 1
1954 (f2cl-lib:array-slice work-%data%
1955 double-float
1956 (n4p5)
1957 ((1
1958 (f2cl-lib:int-add
1959 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
1960 lenaa
1961 nn
1962 1)))
1963 work-%offset%)
1964 1)
1965 (setf pbnprd
1966 (ddot np1
1967 (f2cl-lib:array-slice work-%data%
1968 double-float
1969 (n4p5)
1970 ((1
1971 (f2cl-lib:int-add
1972 (f2cl-lib:int-mul 5
1973 (f2cl-lib:int-add
1974 nn
1975 1))
1976 lenaa
1977 nn
1978 1)))
1979 work-%offset%)
1980 1
1981 (f2cl-lib:array-slice work-%data%
1982 double-float
1983 (n4p5)
1984 ((1
1985 (f2cl-lib:int-add
1986 (f2cl-lib:int-mul 5
1987 (f2cl-lib:int-add
1988 nn
1989 1))
1990 lenaa
1991 nn
1992 1)))
1993 work-%offset%)
1994 1))))
1995 (setf j (f2cl-lib:int-add j 1))
1996 (go label300)
1997 label400
1998 (cond
1999 ((> j imax)
2000 (setf iflag 4)
2001 (go end_label)))
2002 (setf temp
2003 (/
2004 (-
2005 (f2cl-lib:fref work-%data%
2006 ((f2cl-lib:int-add n2p3 nn))
2007 ((1
2008 (f2cl-lib:int-add
2009 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
2010 lenaa
2011 nn
2012 1)))
2013 work-%offset%))
2014 (+ 1.0
2015 (f2cl-lib:fref work-%data%
2016 ((f2cl-lib:int-add np2 nn))
2017 ((1
2018 (f2cl-lib:int-add
2019 (f2cl-lib:int-mul 5
2020 (f2cl-lib:int-add nn 1))
2021 lenaa
2022 nn
2023 1)))
2024 work-%offset%))))
2025 (dcopy np1
2026 (f2cl-lib:array-slice work-%data%
2027 double-float
2028 (n2p3)
2029 ((1
2030 (f2cl-lib:int-add
2031 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
2032 lenaa
2033 nn
2034 1)))
2035 work-%offset%)
2036 1 start 1)
2037 (daxpy np1 temp
2038 (f2cl-lib:array-slice work-%data%
2039 double-float
2040 (np2)
2041 ((1
2042 (f2cl-lib:int-add
2043 (f2cl-lib:int-mul 5 (f2cl-lib:int-add nn 1))
2044 lenaa
2045 nn
2046 1)))
2047 work-%offset%)
2048 1 start 1)
2049 (go end_label)
2050 end_label
2051 (return (values nil nil nil nil nil nil nil nil iflag)))))
2052
2053 (in-package #:cl-user)
2054 #+#.(cl:if (cl:find-package '#:f2cl) '(and) '(or))
2055 (eval-when (:load-toplevel :compile-toplevel :execute)
2056 (setf (gethash 'fortran-to-lisp::pcgns fortran-to-lisp::*f2cl-function-info*)
2057 (fortran-to-lisp::make-f2cl-finfo
2058 :arg-types '((fortran-to-lisp::integer4) (array double-float (*))
2059 (fortran-to-lisp::integer4)
2060 (array fortran-to-lisp::integer4 (*))
2061 (array double-float (*)) (array double-float (*))
2062 (array double-float (*)) (array double-float (*))
2063 (fortran-to-lisp::integer4))
2064 :return-values '(nil nil nil nil nil nil nil nil
2065 fortran-to-lisp::iflag)
2066 :calls '(fortran-to-lisp::ddot fortran-to-lisp::daxpy
2067 fortran-to-lisp::dscal fortran-to-lisp::dnrm2
2068 fortran-to-lisp::dcopy fortran-to-lisp::idamax
2069 fortran-to-lisp::qimuds fortran-to-lisp::multds
2070 fortran-to-lisp::d1mach fortran-to-lisp::mfacds))))
2071
Maxima, a Computer Algebra System