Remove the obsolete DEFMTRFUN-EXTERNAL macro
[maxima.git] / share / hompack / lisp / pcgqs.lisp
blobfd625365a5fd78852fb3a9000046e95985b225d0
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 $")
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))
17 (in-package "HOMPACK")
20 (defun pcgqs (nn aa lenaa maxa pp yp rho start work iflag)
21 (declare (type (array f2cl-lib:integer4 (*)) maxa)
22 (type (array double-float (*)) work start rho yp 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 (yp double-float yp-%data% yp-%offset%)
28 (rho double-float rho-%data% rho-%offset%)
29 (start double-float start-%data% start-%offset%)
30 (work double-float work-%data% work-%offset%)
31 (maxa f2cl-lib:integer4 maxa-%data% maxa-%offset%))
32 (prog ((stillu nil) (stillb nil) (dir 0) (imax 0) (j 0) (lenq 0) (np1 0)
33 (q 0) (res 0) (zb 0) (zu 0) (ab 0.0) (au 0.0) (bb 0.0) (bu 0.0)
34 (dznrm 0.0) (pbnprd 0.0) (punprd 0.0) (rbnprd 0.0) (rbtol 0.0)
35 (rnprd 0.0) (runprd 0.0) (rutol 0.0) (temp 0.0) (zlen 0.0)
36 (ztol 0.0))
37 (declare (type (double-float) ztol zlen temp rutol runprd rnprd rbtol
38 rbnprd punprd pbnprd dznrm bu bb au ab)
39 (type (f2cl-lib:integer4) zu zb res q np1 lenq j imax dir)
40 (type f2cl-lib:logical stillb stillu))
41 (setf np1 (f2cl-lib:int-add nn 1))
42 (setf zu (f2cl-lib:int-add nn 2))
43 (setf zb (f2cl-lib:int-add (f2cl-lib:int-mul 2 nn) 3))
44 (setf res (f2cl-lib:int-add (f2cl-lib:int-mul 3 nn) 4))
45 (setf dir (f2cl-lib:int-add (f2cl-lib:int-mul 4 nn) 5))
46 (setf q (f2cl-lib:int-add (f2cl-lib:int-mul 5 nn) 6))
47 (dcopy lenaa aa 1
48 (f2cl-lib:array-slice work-%data%
49 double-float
50 (q)
51 ((1
52 (f2cl-lib:int-add
53 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
54 lenaa)))
55 work-%offset%)
57 (dcopy np1 yp 1
58 (f2cl-lib:array-slice work-%data%
59 double-float
60 ((+ q lenaa))
61 ((1
62 (f2cl-lib:int-add
63 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
64 lenaa)))
65 work-%offset%)
67 (setf (f2cl-lib:fref maxa-%data%
68 ((f2cl-lib:int-add nn 1))
69 ((1 (f2cl-lib:int-add nn 2)))
70 maxa-%offset%)
71 (f2cl-lib:int-add lenaa 1))
72 (setf (f2cl-lib:fref maxa-%data%
73 ((f2cl-lib:int-add nn 2))
74 ((1 (f2cl-lib:int-add nn 2)))
75 maxa-%offset%)
76 (f2cl-lib:int-add lenaa nn 2))
77 (setf lenq
78 (f2cl-lib:int-sub
79 (f2cl-lib:fref maxa-%data%
80 ((f2cl-lib:int-add nn 2))
81 ((1 (f2cl-lib:int-add nn 2)))
82 maxa-%offset%)
83 1))
84 (gmfads np1
85 (f2cl-lib:array-slice work-%data%
86 double-float
87 (q)
88 ((1
89 (f2cl-lib:int-add
90 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
91 lenaa)))
92 work-%offset%)
93 lenq maxa)
94 (dcopy nn pp 1 work 1)
95 (setf (f2cl-lib:fref work-%data%
96 (np1)
97 ((1
98 (f2cl-lib:int-add
99 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
100 lenaa)))
101 work-%offset%)
102 (coerce 0.0f0 'double-float))
103 (daxpy np1 1.0 yp 1 work 1)
104 (setf imax (f2cl-lib:int-mul 10 np1))
105 (setf stillu f2cl-lib:%true%)
106 (setf stillb f2cl-lib:%true%)
107 (setf ztol (* 100.0f0 (f2cl-lib:d1mach 4)))
108 (setf rutol (* ztol (dnrm2 np1 work 1)))
109 (setf rbtol (* ztol (dnrm2 np1 rho 1)))
110 (multds
111 (f2cl-lib:array-slice work-%data%
112 double-float
113 (res)
115 (f2cl-lib:int-add
116 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
117 lenaa)))
118 work-%offset%)
120 (f2cl-lib:array-slice work-%data%
121 double-float
122 (zu)
124 (f2cl-lib:int-add
125 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
126 lenaa)))
127 work-%offset%)
128 maxa nn lenaa)
129 (setf (f2cl-lib:fref work-%data%
130 ((f2cl-lib:int-add res nn))
132 (f2cl-lib:int-add
133 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
134 lenaa)))
135 work-%offset%)
136 (ddot nn yp 1
137 (f2cl-lib:array-slice work-%data%
138 double-float
139 (zu)
141 (f2cl-lib:int-add
142 (f2cl-lib:int-mul 6
143 (f2cl-lib:int-add nn
145 lenaa)))
146 work-%offset%)
148 (daxpy np1
149 (f2cl-lib:fref work-%data%
150 ((f2cl-lib:int-add zu nn))
152 (f2cl-lib:int-add
153 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
154 lenaa)))
155 work-%offset%)
156 yp 1
157 (f2cl-lib:array-slice work-%data%
158 double-float
159 (res)
161 (f2cl-lib:int-add
162 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
163 lenaa)))
164 work-%offset%)
166 (dscal np1 -1.0
167 (f2cl-lib:array-slice work-%data%
168 double-float
169 (res)
171 (f2cl-lib:int-add
172 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
173 lenaa)))
174 work-%offset%)
176 (daxpy nn -1.0 pp 1
177 (f2cl-lib:array-slice work-%data%
178 double-float
179 (res)
181 (f2cl-lib:int-add
182 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
183 lenaa)))
184 work-%offset%)
186 (daxpy nn -1.0 yp 1
187 (f2cl-lib:array-slice work-%data%
188 double-float
189 (res)
191 (f2cl-lib:int-add
192 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
193 lenaa)))
194 work-%offset%)
196 (solvds np1
197 (f2cl-lib:array-slice work-%data%
198 double-float
201 (f2cl-lib:int-add
202 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
203 lenaa)))
204 work-%offset%)
205 lenq maxa
206 (f2cl-lib:array-slice work-%data%
207 double-float
208 (res)
210 (f2cl-lib:int-add
211 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
212 lenaa)))
213 work-%offset%))
214 (dcopy np1
215 (f2cl-lib:array-slice work-%data%
216 double-float
217 (res)
219 (f2cl-lib:int-add
220 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
221 lenaa)))
222 work-%offset%)
223 1 work 1)
224 (solvds np1
225 (f2cl-lib:array-slice work-%data%
226 double-float
229 (f2cl-lib:int-add
230 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
231 lenaa)))
232 work-%offset%)
233 lenq maxa work)
234 (multds
235 (f2cl-lib:array-slice work-%data%
236 double-float
237 (dir)
239 (f2cl-lib:int-add
240 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
241 lenaa)))
242 work-%offset%)
243 aa work maxa nn lenaa)
244 (setf (f2cl-lib:fref work-%data%
245 ((f2cl-lib:int-add dir nn))
247 (f2cl-lib:int-add
248 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
249 lenaa)))
250 work-%offset%)
251 (ddot nn yp 1 work 1))
252 (daxpy np1
253 (f2cl-lib:fref work-%data%
254 (np1)
256 (f2cl-lib:int-add
257 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
258 lenaa)))
259 work-%offset%)
260 yp 1
261 (f2cl-lib:array-slice work-%data%
262 double-float
263 (dir)
265 (f2cl-lib:int-add
266 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
267 lenaa)))
268 work-%offset%)
270 (setf runprd
271 (ddot np1
272 (f2cl-lib:array-slice work-%data%
273 double-float
274 (res)
276 (f2cl-lib:int-add
277 (f2cl-lib:int-mul 6
278 (f2cl-lib:int-add nn
280 lenaa)))
281 work-%offset%)
283 (f2cl-lib:array-slice work-%data%
284 double-float
285 (res)
287 (f2cl-lib:int-add
288 (f2cl-lib:int-mul 6
289 (f2cl-lib:int-add nn
291 lenaa)))
292 work-%offset%)
294 (setf punprd
295 (ddot np1
296 (f2cl-lib:array-slice work-%data%
297 double-float
298 (dir)
300 (f2cl-lib:int-add
301 (f2cl-lib:int-mul 6
302 (f2cl-lib:int-add nn
304 lenaa)))
305 work-%offset%)
307 (f2cl-lib:array-slice work-%data%
308 double-float
309 (dir)
311 (f2cl-lib:int-add
312 (f2cl-lib:int-mul 6
313 (f2cl-lib:int-add nn
315 lenaa)))
316 work-%offset%)
318 (setf j 1)
319 label100
320 (if (not (and stillu (<= j imax))) (go label200))
321 (cond
322 ((> (f2cl-lib:fsqrt runprd) rutol)
323 (cond
324 ((= punprd 0.0f0)
325 (multds
326 (f2cl-lib:array-slice work-%data%
327 double-float
328 (res)
330 (f2cl-lib:int-add
331 (f2cl-lib:int-mul 6
332 (f2cl-lib:int-add nn
334 lenaa)))
335 work-%offset%)
337 (f2cl-lib:array-slice work-%data%
338 double-float
339 (zu)
341 (f2cl-lib:int-add
342 (f2cl-lib:int-mul 6
343 (f2cl-lib:int-add nn
345 lenaa)))
346 work-%offset%)
347 maxa nn lenaa)
348 (setf (f2cl-lib:fref work-%data%
349 ((f2cl-lib:int-add res nn))
351 (f2cl-lib:int-add
352 (f2cl-lib:int-mul 6
353 (f2cl-lib:int-add nn 1))
354 lenaa)))
355 work-%offset%)
356 (ddot nn yp 1
357 (f2cl-lib:array-slice work-%data%
358 double-float
359 (zu)
361 (f2cl-lib:int-add
362 (f2cl-lib:int-mul 6
363 (f2cl-lib:int-add
366 lenaa)))
367 work-%offset%)
369 (daxpy np1
370 (f2cl-lib:fref work-%data%
371 ((f2cl-lib:int-add zu nn))
373 (f2cl-lib:int-add
374 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
375 lenaa)))
376 work-%offset%)
377 yp 1
378 (f2cl-lib:array-slice work-%data%
379 double-float
380 (res)
382 (f2cl-lib:int-add
383 (f2cl-lib:int-mul 6
384 (f2cl-lib:int-add nn
386 lenaa)))
387 work-%offset%)
389 (dscal np1 -1.0
390 (f2cl-lib:array-slice work-%data%
391 double-float
392 (res)
394 (f2cl-lib:int-add
395 (f2cl-lib:int-mul 6
396 (f2cl-lib:int-add nn
398 lenaa)))
399 work-%offset%)
401 (daxpy nn -1.0 pp 1
402 (f2cl-lib:array-slice work-%data%
403 double-float
404 (res)
406 (f2cl-lib:int-add
407 (f2cl-lib:int-mul 6
408 (f2cl-lib:int-add nn
410 lenaa)))
411 work-%offset%)
413 (daxpy nn -1.0 yp 1
414 (f2cl-lib:array-slice work-%data%
415 double-float
416 (res)
418 (f2cl-lib:int-add
419 (f2cl-lib:int-mul 6
420 (f2cl-lib:int-add nn
422 lenaa)))
423 work-%offset%)
425 (solvds np1
426 (f2cl-lib:array-slice work-%data%
427 double-float
430 (f2cl-lib:int-add
431 (f2cl-lib:int-mul 6
432 (f2cl-lib:int-add nn
434 lenaa)))
435 work-%offset%)
436 lenq maxa
437 (f2cl-lib:array-slice work-%data%
438 double-float
439 (res)
441 (f2cl-lib:int-add
442 (f2cl-lib:int-mul 6
443 (f2cl-lib:int-add nn
445 lenaa)))
446 work-%offset%))
447 (dcopy np1
448 (f2cl-lib:array-slice work-%data%
449 double-float
450 (res)
452 (f2cl-lib:int-add
453 (f2cl-lib:int-mul 6
454 (f2cl-lib:int-add nn
456 lenaa)))
457 work-%offset%)
458 1 work 1)
459 (solvds np1
460 (f2cl-lib:array-slice work-%data%
461 double-float
464 (f2cl-lib:int-add
465 (f2cl-lib:int-mul 6
466 (f2cl-lib:int-add nn
468 lenaa)))
469 work-%offset%)
470 lenq maxa work)
471 (multds
472 (f2cl-lib:array-slice work-%data%
473 double-float
474 (dir)
476 (f2cl-lib:int-add
477 (f2cl-lib:int-mul 6
478 (f2cl-lib:int-add nn
480 lenaa)))
481 work-%offset%)
482 aa work maxa nn lenaa)
483 (setf (f2cl-lib:fref work-%data%
484 ((f2cl-lib:int-add dir nn))
486 (f2cl-lib:int-add
487 (f2cl-lib:int-mul 6
488 (f2cl-lib:int-add nn 1))
489 lenaa)))
490 work-%offset%)
491 (ddot nn yp 1 work 1))
492 (daxpy np1
493 (f2cl-lib:fref work-%data%
494 (np1)
496 (f2cl-lib:int-add
497 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
498 lenaa)))
499 work-%offset%)
500 yp 1
501 (f2cl-lib:array-slice work-%data%
502 double-float
503 (dir)
505 (f2cl-lib:int-add
506 (f2cl-lib:int-mul 6
507 (f2cl-lib:int-add nn
509 lenaa)))
510 work-%offset%)
512 (setf runprd
513 (ddot np1
514 (f2cl-lib:array-slice work-%data%
515 double-float
516 (res)
518 (f2cl-lib:int-add
519 (f2cl-lib:int-mul 6
520 (f2cl-lib:int-add
523 lenaa)))
524 work-%offset%)
526 (f2cl-lib:array-slice work-%data%
527 double-float
528 (res)
530 (f2cl-lib:int-add
531 (f2cl-lib:int-mul 6
532 (f2cl-lib:int-add
535 lenaa)))
536 work-%offset%)
538 (setf punprd
539 (ddot np1
540 (f2cl-lib:array-slice work-%data%
541 double-float
542 (dir)
544 (f2cl-lib:int-add
545 (f2cl-lib:int-mul 6
546 (f2cl-lib:int-add
549 lenaa)))
550 work-%offset%)
552 (f2cl-lib:array-slice work-%data%
553 double-float
554 (dir)
556 (f2cl-lib:int-add
557 (f2cl-lib:int-mul 6
558 (f2cl-lib:int-add
561 lenaa)))
562 work-%offset%)
564 (cond
565 ((<= (f2cl-lib:fsqrt runprd) rutol)
566 (setf stillu f2cl-lib:%false%)))))
567 (cond
568 (stillu
569 (setf au (/ runprd punprd))
570 (daxpy np1 au
571 (f2cl-lib:array-slice work-%data%
572 double-float
573 (dir)
575 (f2cl-lib:int-add
576 (f2cl-lib:int-mul 6
577 (f2cl-lib:int-add nn
579 lenaa)))
580 work-%offset%)
582 (f2cl-lib:array-slice work-%data%
583 double-float
584 (zu)
586 (f2cl-lib:int-add
587 (f2cl-lib:int-mul 6
588 (f2cl-lib:int-add nn
590 lenaa)))
591 work-%offset%)
593 (setf dznrm (* au (f2cl-lib:fsqrt punprd)))
594 (setf zlen
595 (dnrm2 np1
596 (f2cl-lib:array-slice work-%data%
597 double-float
598 (zu)
600 (f2cl-lib:int-add
601 (f2cl-lib:int-mul 6
602 (f2cl-lib:int-add
605 lenaa)))
606 work-%offset%)
608 (if (< (/ dznrm zlen) ztol) (setf stillu f2cl-lib:%false%)))))
610 (setf stillu f2cl-lib:%false%)))
611 (cond
612 (stillu
613 (multds work aa
614 (f2cl-lib:array-slice work-%data%
615 double-float
616 (dir)
618 (f2cl-lib:int-add
619 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
620 lenaa)))
621 work-%offset%)
622 maxa nn lenaa)
623 (setf (f2cl-lib:fref work-%data%
624 (np1)
626 (f2cl-lib:int-add
627 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
628 lenaa)))
629 work-%offset%)
630 (ddot nn yp 1
631 (f2cl-lib:array-slice work-%data%
632 double-float
633 (dir)
635 (f2cl-lib:int-add
636 (f2cl-lib:int-mul 6
637 (f2cl-lib:int-add
640 lenaa)))
641 work-%offset%)
643 (daxpy np1
644 (f2cl-lib:fref work-%data%
645 ((f2cl-lib:int-add dir nn))
647 (f2cl-lib:int-add
648 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
649 lenaa)))
650 work-%offset%)
651 yp 1 work 1)
652 (solvds np1
653 (f2cl-lib:array-slice work-%data%
654 double-float
657 (f2cl-lib:int-add
658 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
659 lenaa)))
660 work-%offset%)
661 lenq maxa work)
662 (daxpy np1 (- au) work 1
663 (f2cl-lib:array-slice work-%data%
664 double-float
665 (res)
667 (f2cl-lib:int-add
668 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
669 lenaa)))
670 work-%offset%)
672 (setf rnprd
673 (ddot np1
674 (f2cl-lib:array-slice work-%data%
675 double-float
676 (res)
678 (f2cl-lib:int-add
679 (f2cl-lib:int-mul 6
680 (f2cl-lib:int-add
683 lenaa)))
684 work-%offset%)
686 (f2cl-lib:array-slice work-%data%
687 double-float
688 (res)
690 (f2cl-lib:int-add
691 (f2cl-lib:int-mul 6
692 (f2cl-lib:int-add
695 lenaa)))
696 work-%offset%)
698 (setf bu (/ rnprd runprd))
699 (setf runprd rnprd)
700 (dcopy np1
701 (f2cl-lib:array-slice work-%data%
702 double-float
703 (res)
705 (f2cl-lib:int-add
706 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
707 lenaa)))
708 work-%offset%)
709 1 work 1)
710 (solvds np1
711 (f2cl-lib:array-slice work-%data%
712 double-float
715 (f2cl-lib:int-add
716 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
717 lenaa)))
718 work-%offset%)
719 lenq maxa work)
720 (multds start aa work maxa nn lenaa)
721 (setf (f2cl-lib:fref start-%data%
722 (np1)
723 ((1 (f2cl-lib:int-add nn 1)))
724 start-%offset%)
725 (ddot nn yp 1 work 1))
726 (daxpy np1
727 (f2cl-lib:fref work-%data%
728 (np1)
730 (f2cl-lib:int-add
731 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
732 lenaa)))
733 work-%offset%)
734 yp 1 start 1)
735 (daxpy np1 bu
736 (f2cl-lib:array-slice work-%data%
737 double-float
738 (dir)
740 (f2cl-lib:int-add
741 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
742 lenaa)))
743 work-%offset%)
744 1 start 1)
745 (dcopy np1 start 1
746 (f2cl-lib:array-slice work-%data%
747 double-float
748 (dir)
750 (f2cl-lib:int-add
751 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
752 lenaa)))
753 work-%offset%)
755 (setf punprd
756 (ddot np1
757 (f2cl-lib:array-slice work-%data%
758 double-float
759 (dir)
761 (f2cl-lib:int-add
762 (f2cl-lib:int-mul 6
763 (f2cl-lib:int-add
766 lenaa)))
767 work-%offset%)
769 (f2cl-lib:array-slice work-%data%
770 double-float
771 (dir)
773 (f2cl-lib:int-add
774 (f2cl-lib:int-mul 6
775 (f2cl-lib:int-add
778 lenaa)))
779 work-%offset%)
780 1))))
781 (setf j (f2cl-lib:int-add j 1))
782 (go label100)
783 label200
784 (cond
785 ((> j imax)
786 (setf iflag 4)
787 (go end_label)))
788 (multds
789 (f2cl-lib:array-slice work-%data%
790 double-float
791 (res)
793 (f2cl-lib:int-add
794 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
795 lenaa)))
796 work-%offset%)
798 (f2cl-lib:array-slice work-%data%
799 double-float
800 (zb)
802 (f2cl-lib:int-add
803 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
804 lenaa)))
805 work-%offset%)
806 maxa nn lenaa)
807 (setf (f2cl-lib:fref work-%data%
808 ((f2cl-lib:int-add res nn))
810 (f2cl-lib:int-add
811 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
812 lenaa)))
813 work-%offset%)
814 (ddot nn yp 1
815 (f2cl-lib:array-slice work-%data%
816 double-float
817 (zb)
819 (f2cl-lib:int-add
820 (f2cl-lib:int-mul 6
821 (f2cl-lib:int-add nn
823 lenaa)))
824 work-%offset%)
826 (daxpy np1
827 (f2cl-lib:fref work-%data%
828 ((f2cl-lib:int-add zb nn))
830 (f2cl-lib:int-add
831 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
832 lenaa)))
833 work-%offset%)
834 yp 1
835 (f2cl-lib:array-slice work-%data%
836 double-float
837 (res)
839 (f2cl-lib:int-add
840 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
841 lenaa)))
842 work-%offset%)
844 (dscal np1 -1.0
845 (f2cl-lib:array-slice work-%data%
846 double-float
847 (res)
849 (f2cl-lib:int-add
850 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
851 lenaa)))
852 work-%offset%)
854 (daxpy np1 -1.0 rho 1
855 (f2cl-lib:array-slice work-%data%
856 double-float
857 (res)
859 (f2cl-lib:int-add
860 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
861 lenaa)))
862 work-%offset%)
864 (solvds np1
865 (f2cl-lib:array-slice work-%data%
866 double-float
869 (f2cl-lib:int-add
870 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
871 lenaa)))
872 work-%offset%)
873 lenq maxa
874 (f2cl-lib:array-slice work-%data%
875 double-float
876 (res)
878 (f2cl-lib:int-add
879 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
880 lenaa)))
881 work-%offset%))
882 (dcopy np1
883 (f2cl-lib:array-slice work-%data%
884 double-float
885 (res)
887 (f2cl-lib:int-add
888 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
889 lenaa)))
890 work-%offset%)
891 1 work 1)
892 (solvds np1
893 (f2cl-lib:array-slice work-%data%
894 double-float
897 (f2cl-lib:int-add
898 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
899 lenaa)))
900 work-%offset%)
901 lenq maxa work)
902 (multds
903 (f2cl-lib:array-slice work-%data%
904 double-float
905 (dir)
907 (f2cl-lib:int-add
908 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
909 lenaa)))
910 work-%offset%)
911 aa work maxa nn lenaa)
912 (setf (f2cl-lib:fref work-%data%
913 ((f2cl-lib:int-add dir nn))
915 (f2cl-lib:int-add
916 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
917 lenaa)))
918 work-%offset%)
919 (ddot nn yp 1 work 1))
920 (daxpy np1
921 (f2cl-lib:fref work-%data%
922 (np1)
924 (f2cl-lib:int-add
925 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
926 lenaa)))
927 work-%offset%)
928 yp 1
929 (f2cl-lib:array-slice work-%data%
930 double-float
931 (dir)
933 (f2cl-lib:int-add
934 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
935 lenaa)))
936 work-%offset%)
938 (setf rbnprd
939 (ddot np1
940 (f2cl-lib:array-slice work-%data%
941 double-float
942 (res)
944 (f2cl-lib:int-add
945 (f2cl-lib:int-mul 6
946 (f2cl-lib:int-add nn
948 lenaa)))
949 work-%offset%)
951 (f2cl-lib:array-slice work-%data%
952 double-float
953 (res)
955 (f2cl-lib:int-add
956 (f2cl-lib:int-mul 6
957 (f2cl-lib:int-add nn
959 lenaa)))
960 work-%offset%)
962 (setf pbnprd
963 (ddot np1
964 (f2cl-lib:array-slice work-%data%
965 double-float
966 (dir)
968 (f2cl-lib:int-add
969 (f2cl-lib:int-mul 6
970 (f2cl-lib:int-add nn
972 lenaa)))
973 work-%offset%)
975 (f2cl-lib:array-slice work-%data%
976 double-float
977 (dir)
979 (f2cl-lib:int-add
980 (f2cl-lib:int-mul 6
981 (f2cl-lib:int-add nn
983 lenaa)))
984 work-%offset%)
986 (setf j 1)
987 label300
988 (if (not (and stillb (<= j imax))) (go label400))
989 (cond
990 ((> (f2cl-lib:fsqrt rbnprd) rbtol)
991 (cond
992 ((= pbnprd 0.0f0)
993 (multds
994 (f2cl-lib:array-slice work-%data%
995 double-float
996 (res)
998 (f2cl-lib:int-add
999 (f2cl-lib:int-mul 6
1000 (f2cl-lib:int-add nn
1002 lenaa)))
1003 work-%offset%)
1005 (f2cl-lib:array-slice work-%data%
1006 double-float
1007 (zb)
1009 (f2cl-lib:int-add
1010 (f2cl-lib:int-mul 6
1011 (f2cl-lib:int-add nn
1013 lenaa)))
1014 work-%offset%)
1015 maxa nn lenaa)
1016 (setf (f2cl-lib:fref work-%data%
1017 ((f2cl-lib:int-add res nn))
1019 (f2cl-lib:int-add
1020 (f2cl-lib:int-mul 6
1021 (f2cl-lib:int-add nn 1))
1022 lenaa)))
1023 work-%offset%)
1024 (ddot nn yp 1
1025 (f2cl-lib:array-slice work-%data%
1026 double-float
1027 (zb)
1029 (f2cl-lib:int-add
1030 (f2cl-lib:int-mul 6
1031 (f2cl-lib:int-add
1034 lenaa)))
1035 work-%offset%)
1037 (daxpy np1
1038 (f2cl-lib:fref work-%data%
1039 ((f2cl-lib:int-add zb nn))
1041 (f2cl-lib:int-add
1042 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
1043 lenaa)))
1044 work-%offset%)
1045 yp 1
1046 (f2cl-lib:array-slice work-%data%
1047 double-float
1048 (res)
1050 (f2cl-lib:int-add
1051 (f2cl-lib:int-mul 6
1052 (f2cl-lib:int-add nn
1054 lenaa)))
1055 work-%offset%)
1057 (dscal np1 -1.0
1058 (f2cl-lib:array-slice work-%data%
1059 double-float
1060 (res)
1062 (f2cl-lib:int-add
1063 (f2cl-lib:int-mul 6
1064 (f2cl-lib:int-add nn
1066 lenaa)))
1067 work-%offset%)
1069 (daxpy np1 -1.0 rho 1
1070 (f2cl-lib:array-slice work-%data%
1071 double-float
1072 (res)
1074 (f2cl-lib:int-add
1075 (f2cl-lib:int-mul 6
1076 (f2cl-lib:int-add nn
1078 lenaa)))
1079 work-%offset%)
1081 (solvds np1
1082 (f2cl-lib:array-slice work-%data%
1083 double-float
1086 (f2cl-lib:int-add
1087 (f2cl-lib:int-mul 6
1088 (f2cl-lib:int-add nn
1090 lenaa)))
1091 work-%offset%)
1092 lenq maxa
1093 (f2cl-lib:array-slice work-%data%
1094 double-float
1095 (res)
1097 (f2cl-lib:int-add
1098 (f2cl-lib:int-mul 6
1099 (f2cl-lib:int-add nn
1101 lenaa)))
1102 work-%offset%))
1103 (dcopy np1
1104 (f2cl-lib:array-slice work-%data%
1105 double-float
1106 (res)
1108 (f2cl-lib:int-add
1109 (f2cl-lib:int-mul 6
1110 (f2cl-lib:int-add nn
1112 lenaa)))
1113 work-%offset%)
1114 1 work 1)
1115 (solvds np1
1116 (f2cl-lib:array-slice work-%data%
1117 double-float
1120 (f2cl-lib:int-add
1121 (f2cl-lib:int-mul 6
1122 (f2cl-lib:int-add nn
1124 lenaa)))
1125 work-%offset%)
1126 lenq maxa work)
1127 (multds
1128 (f2cl-lib:array-slice work-%data%
1129 double-float
1130 (dir)
1132 (f2cl-lib:int-add
1133 (f2cl-lib:int-mul 6
1134 (f2cl-lib:int-add nn
1136 lenaa)))
1137 work-%offset%)
1138 aa work maxa nn lenaa)
1139 (setf (f2cl-lib:fref work-%data%
1140 ((f2cl-lib:int-add dir nn))
1142 (f2cl-lib:int-add
1143 (f2cl-lib:int-mul 6
1144 (f2cl-lib:int-add nn 1))
1145 lenaa)))
1146 work-%offset%)
1147 (ddot nn yp 1 work 1))
1148 (daxpy np1
1149 (f2cl-lib:fref work-%data%
1150 (np1)
1152 (f2cl-lib:int-add
1153 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
1154 lenaa)))
1155 work-%offset%)
1156 yp 1
1157 (f2cl-lib:array-slice work-%data%
1158 double-float
1159 (dir)
1161 (f2cl-lib:int-add
1162 (f2cl-lib:int-mul 6
1163 (f2cl-lib:int-add nn
1165 lenaa)))
1166 work-%offset%)
1168 (setf rbnprd
1169 (ddot np1
1170 (f2cl-lib:array-slice work-%data%
1171 double-float
1172 (res)
1174 (f2cl-lib:int-add
1175 (f2cl-lib:int-mul 6
1176 (f2cl-lib:int-add
1179 lenaa)))
1180 work-%offset%)
1182 (f2cl-lib:array-slice work-%data%
1183 double-float
1184 (res)
1186 (f2cl-lib:int-add
1187 (f2cl-lib:int-mul 6
1188 (f2cl-lib:int-add
1191 lenaa)))
1192 work-%offset%)
1194 (setf pbnprd
1195 (ddot np1
1196 (f2cl-lib:array-slice work-%data%
1197 double-float
1198 (dir)
1200 (f2cl-lib:int-add
1201 (f2cl-lib:int-mul 6
1202 (f2cl-lib:int-add
1205 lenaa)))
1206 work-%offset%)
1208 (f2cl-lib:array-slice work-%data%
1209 double-float
1210 (dir)
1212 (f2cl-lib:int-add
1213 (f2cl-lib:int-mul 6
1214 (f2cl-lib:int-add
1217 lenaa)))
1218 work-%offset%)
1220 (cond
1221 ((<= (f2cl-lib:fsqrt rbnprd) rbtol)
1222 (setf stillb f2cl-lib:%false%)))))
1223 (cond
1224 (stillb
1225 (setf ab (/ rbnprd pbnprd))
1226 (daxpy np1 ab
1227 (f2cl-lib:array-slice work-%data%
1228 double-float
1229 (dir)
1231 (f2cl-lib:int-add
1232 (f2cl-lib:int-mul 6
1233 (f2cl-lib:int-add nn
1235 lenaa)))
1236 work-%offset%)
1238 (f2cl-lib:array-slice work-%data%
1239 double-float
1240 (zb)
1242 (f2cl-lib:int-add
1243 (f2cl-lib:int-mul 6
1244 (f2cl-lib:int-add nn
1246 lenaa)))
1247 work-%offset%)
1249 (setf dznrm (* ab (f2cl-lib:fsqrt pbnprd)))
1250 (setf zlen
1251 (dnrm2 np1
1252 (f2cl-lib:array-slice work-%data%
1253 double-float
1254 (zb)
1256 (f2cl-lib:int-add
1257 (f2cl-lib:int-mul 6
1258 (f2cl-lib:int-add
1261 lenaa)))
1262 work-%offset%)
1264 (if (< (/ dznrm zlen) ztol) (setf stillb f2cl-lib:%false%)))))
1266 (setf stillb f2cl-lib:%false%)))
1267 (cond
1268 (stillb
1269 (multds work aa
1270 (f2cl-lib:array-slice work-%data%
1271 double-float
1272 (dir)
1274 (f2cl-lib:int-add
1275 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
1276 lenaa)))
1277 work-%offset%)
1278 maxa nn lenaa)
1279 (setf (f2cl-lib:fref work-%data%
1280 (np1)
1282 (f2cl-lib:int-add
1283 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
1284 lenaa)))
1285 work-%offset%)
1286 (ddot nn yp 1
1287 (f2cl-lib:array-slice work-%data%
1288 double-float
1289 (dir)
1291 (f2cl-lib:int-add
1292 (f2cl-lib:int-mul 6
1293 (f2cl-lib:int-add
1296 lenaa)))
1297 work-%offset%)
1299 (daxpy np1
1300 (f2cl-lib:fref work-%data%
1301 ((f2cl-lib:int-add dir nn))
1303 (f2cl-lib:int-add
1304 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
1305 lenaa)))
1306 work-%offset%)
1307 yp 1 work 1)
1308 (solvds np1
1309 (f2cl-lib:array-slice work-%data%
1310 double-float
1313 (f2cl-lib:int-add
1314 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
1315 lenaa)))
1316 work-%offset%)
1317 lenq maxa work)
1318 (daxpy np1 (- ab) work 1
1319 (f2cl-lib:array-slice work-%data%
1320 double-float
1321 (res)
1323 (f2cl-lib:int-add
1324 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
1325 lenaa)))
1326 work-%offset%)
1328 (setf rnprd
1329 (ddot np1
1330 (f2cl-lib:array-slice work-%data%
1331 double-float
1332 (res)
1334 (f2cl-lib:int-add
1335 (f2cl-lib:int-mul 6
1336 (f2cl-lib:int-add
1339 lenaa)))
1340 work-%offset%)
1342 (f2cl-lib:array-slice work-%data%
1343 double-float
1344 (res)
1346 (f2cl-lib:int-add
1347 (f2cl-lib:int-mul 6
1348 (f2cl-lib:int-add
1351 lenaa)))
1352 work-%offset%)
1354 (setf bb (/ rnprd rbnprd))
1355 (setf rbnprd rnprd)
1356 (dcopy np1
1357 (f2cl-lib:array-slice work-%data%
1358 double-float
1359 (res)
1361 (f2cl-lib:int-add
1362 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
1363 lenaa)))
1364 work-%offset%)
1365 1 work 1)
1366 (solvds np1
1367 (f2cl-lib:array-slice work-%data%
1368 double-float
1371 (f2cl-lib:int-add
1372 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
1373 lenaa)))
1374 work-%offset%)
1375 lenq maxa work)
1376 (multds start aa work maxa nn lenaa)
1377 (setf (f2cl-lib:fref start-%data%
1378 (np1)
1379 ((1 (f2cl-lib:int-add nn 1)))
1380 start-%offset%)
1381 (ddot nn yp 1 work 1))
1382 (daxpy np1
1383 (f2cl-lib:fref work-%data%
1384 (np1)
1386 (f2cl-lib:int-add
1387 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
1388 lenaa)))
1389 work-%offset%)
1390 yp 1 start 1)
1391 (daxpy np1 bb
1392 (f2cl-lib:array-slice work-%data%
1393 double-float
1394 (dir)
1396 (f2cl-lib:int-add
1397 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
1398 lenaa)))
1399 work-%offset%)
1400 1 start 1)
1401 (dcopy np1 start 1
1402 (f2cl-lib:array-slice work-%data%
1403 double-float
1404 (dir)
1406 (f2cl-lib:int-add
1407 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
1408 lenaa)))
1409 work-%offset%)
1411 (setf pbnprd
1412 (ddot np1
1413 (f2cl-lib:array-slice work-%data%
1414 double-float
1415 (dir)
1417 (f2cl-lib:int-add
1418 (f2cl-lib:int-mul 6
1419 (f2cl-lib:int-add
1422 lenaa)))
1423 work-%offset%)
1425 (f2cl-lib:array-slice work-%data%
1426 double-float
1427 (dir)
1429 (f2cl-lib:int-add
1430 (f2cl-lib:int-mul 6
1431 (f2cl-lib:int-add
1434 lenaa)))
1435 work-%offset%)
1436 1))))
1437 (setf j (f2cl-lib:int-add j 1))
1438 (go label300)
1439 label400
1440 (cond
1441 ((> j imax)
1442 (setf iflag 4)
1443 (go end_label)))
1444 (setf temp
1447 (f2cl-lib:fref work-%data%
1448 ((f2cl-lib:int-add zb nn))
1450 (f2cl-lib:int-add
1451 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
1452 lenaa)))
1453 work-%offset%))
1454 (+ 1.0
1455 (f2cl-lib:fref work-%data%
1456 ((f2cl-lib:int-add zu nn))
1458 (f2cl-lib:int-add
1459 (f2cl-lib:int-mul 6
1460 (f2cl-lib:int-add nn 1))
1461 lenaa)))
1462 work-%offset%))))
1463 (dcopy np1
1464 (f2cl-lib:array-slice work-%data%
1465 double-float
1466 (zb)
1468 (f2cl-lib:int-add
1469 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
1470 lenaa)))
1471 work-%offset%)
1472 1 start 1)
1473 (daxpy np1 temp
1474 (f2cl-lib:array-slice work-%data%
1475 double-float
1476 (zu)
1478 (f2cl-lib:int-add
1479 (f2cl-lib:int-mul 6 (f2cl-lib:int-add nn 1))
1480 lenaa)))
1481 work-%offset%)
1482 1 start 1)
1483 (go end_label)
1484 end_label
1485 (return (values nil nil nil nil nil nil nil nil nil iflag)))))
1487 (in-package #:cl-user)
1488 #+#.(cl:if (cl:find-package '#:f2cl) '(and) '(or))
1489 (eval-when (:load-toplevel :compile-toplevel :execute)
1490 (setf (gethash 'fortran-to-lisp::pcgqs fortran-to-lisp::*f2cl-function-info*)
1491 (fortran-to-lisp::make-f2cl-finfo
1492 :arg-types '((fortran-to-lisp::integer4) (array double-float (*))
1493 (fortran-to-lisp::integer4)
1494 (array fortran-to-lisp::integer4 (*))
1495 (array double-float (*)) (array double-float (*))
1496 (array double-float (*)) (array double-float (*))
1497 (array double-float (*)) (fortran-to-lisp::integer4))
1498 :return-values '(nil nil nil nil nil nil nil nil nil
1499 fortran-to-lisp::iflag)
1500 :calls '(fortran-to-lisp::dscal fortran-to-lisp::ddot
1501 fortran-to-lisp::dnrm2 fortran-to-lisp::daxpy
1502 fortran-to-lisp::dcopy fortran-to-lisp::solvds
1503 fortran-to-lisp::multds fortran-to-lisp::d1mach
1504 fortran-to-lisp::gmfads))))