MEVALP_TR: return result of MEVALP1_TR instead of unknown
[maxima.git] / share / hompack / lisp / rootns.lisp
bloba6d67a3d1a5d392f86d3d68b5c6a1730ac414991
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 (let* ((limit 20))
21 (declare (type (f2cl-lib:integer4 20 20) limit) (ignorable limit))
22 (defun rootns
23 (n nfe iflag relerr abserr y yp yold ypold a qr lenqr pivot work par
24 ipar)
25 (declare (type (array f2cl-lib:integer4 (*)) ipar)
26 (type (array double-float (*)) par)
27 (type (array f2cl-lib:integer4 (*)) pivot)
28 (type (array double-float (*)) work qr a ypold yold yp y)
29 (type (double-float) abserr relerr)
30 (type (f2cl-lib:integer4) lenqr iflag nfe n))
31 (f2cl-lib:with-multi-array-data
32 ((y double-float y-%data% y-%offset%)
33 (yp double-float yp-%data% yp-%offset%)
34 (yold double-float yold-%data% yold-%offset%)
35 (ypold double-float ypold-%data% ypold-%offset%)
36 (a double-float a-%data% a-%offset%)
37 (qr double-float qr-%data% qr-%offset%)
38 (work double-float work-%data% work-%offset%)
39 (pivot f2cl-lib:integer4 pivot-%data% pivot-%offset%)
40 (par double-float par-%data% par-%offset%)
41 (ipar f2cl-lib:integer4 ipar-%data% ipar-%offset%))
42 (labels ((dd01 (f0 f1 dels)
43 (f2cl-lib:f2cl/ (+ f1 (- f0)) dels))
44 (dd001 (f0 fp0 f1 dels)
45 (/ (- (dd01 f0 f1 dels) fp0) dels))
46 (dd011 (f0 f1 fp1 dels)
47 (/ (- fp1 (dd01 f0 f1 dels)) dels))
48 (dd0011 (f0 fp0 f1 fp1 dels)
49 (/ (- (dd011 f0 f1 fp1 dels) (dd001 f0 fp0 f1 dels)) dels))
50 (qofs (f0 fp0 f1 fp1 dels s)
55 (+ (* (dd0011 f0 fp0 f1 fp1 dels) (- s dels))
56 (dd001 f0 fp0 f1 dels))
58 fp0)
60 f0)))
61 (declare (ftype (function (double-float double-float double-float)
62 (values double-float &rest t))
63 dd01))
64 (declare (ftype (function
65 (double-float double-float double-float double-float)
66 (values double-float &rest t))
67 dd001))
68 (declare (ftype (function
69 (double-float double-float double-float double-float)
70 (values double-float &rest t))
71 dd011))
72 (declare (ftype (function
73 (double-float double-float double-float double-float
74 double-float)
75 (values double-float &rest t))
76 dd0011))
77 (declare (ftype (function
78 (double-float double-float double-float double-float
79 double-float double-float)
80 (values double-float &rest t))
81 qofs))
82 (prog ((ipp 0) (irho 0) (itangw 0) (itz 0) (iw 0) (iwp 0) (iz0 0)
83 (iz1 0) (judy 0) (jw 0) (lcode 0) (np1 0) (aerr 0.0) (dels 0.0)
84 (f0 0.0) (f1 0.0) (fp0 0.0) (fp1 0.0) (qsout 0.0) (rerr 0.0)
85 (s 0.0) (sa 0.0) (sb 0.0) (sout 0.0) (u 0.0))
86 (declare (type (f2cl-lib:integer4) ipp irho itangw itz iw iwp iz0 iz1
87 judy jw lcode np1)
88 (type (double-float) aerr dels f0 f1 fp0 fp1 qsout rerr s sa
89 sb sout u))
90 (setf u (f2cl-lib:d1mach 4))
91 (setf rerr (max relerr u))
92 (setf aerr (max abserr 0.0))
93 (setf np1 (f2cl-lib:int-add n 1))
94 (setf ipp 1)
95 (setf irho (f2cl-lib:int-add n 1))
96 (setf iw (f2cl-lib:int-add irho n))
97 (setf iwp (f2cl-lib:int-add iw np1))
98 (setf itz (f2cl-lib:int-add iwp np1))
99 (setf iz0 (f2cl-lib:int-add itz np1))
100 (setf iz1 (f2cl-lib:int-add iz0 np1))
101 (setf itangw (f2cl-lib:int-add iz1 np1))
102 label100
103 (f2cl-lib:fdo (judy 1 (f2cl-lib:int-add judy 1))
104 ((> judy limit) nil)
105 (tagbody
106 (f2cl-lib:fdo (jw 1 (f2cl-lib:int-add jw 1))
107 ((> jw np1) nil)
108 (tagbody
109 (setf (f2cl-lib:fref work-%data%
110 ((f2cl-lib:int-sub
111 (f2cl-lib:int-add itz jw)
114 (f2cl-lib:int-add
115 (f2cl-lib:int-mul 13
116 (f2cl-lib:int-add n
118 (f2cl-lib:int-mul 2 n)
119 lenqr)))
120 work-%offset%)
122 (f2cl-lib:fref y-%data%
123 (jw)
124 ((1 (f2cl-lib:int-add n 1)))
125 y-%offset%)
126 (f2cl-lib:fref yold-%data%
127 (jw)
128 ((1 (f2cl-lib:int-add n 1)))
129 yold-%offset%)))
130 label110))
131 (setf dels
132 (dnrm2 np1
133 (f2cl-lib:array-slice work-%data%
134 double-float
135 (itz)
137 (f2cl-lib:int-add
138 (f2cl-lib:int-mul 13
139 (f2cl-lib:int-add
142 (f2cl-lib:int-mul 2 n)
143 lenqr)))
144 work-%offset%)
146 (setf sa (coerce 0.0f0 'double-float))
147 (setf sb dels)
148 (setf lcode 1)
149 label130
150 (multiple-value-bind (var-0 var-1 var-2 var-3 var-4 var-5 var-6)
151 (root sout qsout sa sb rerr aerr lcode)
152 (declare (ignore var-1 var-4 var-5))
153 (setf sout var-0)
154 (setf sa var-2)
155 (setf sb var-3)
156 (setf lcode var-6))
157 (if (> lcode 0) (go label140))
158 (setf qsout
160 (qofs
161 (f2cl-lib:fref yold-%data%
162 (np1)
163 ((1 (f2cl-lib:int-add n 1)))
164 yold-%offset%)
165 (f2cl-lib:fref ypold-%data%
166 (np1)
167 ((1 (f2cl-lib:int-add n 1)))
168 ypold-%offset%)
169 (f2cl-lib:fref y-%data%
170 (np1)
171 ((1 (f2cl-lib:int-add n 1)))
172 y-%offset%)
173 (f2cl-lib:fref yp-%data%
174 (np1)
175 ((1 (f2cl-lib:int-add n 1)))
176 yp-%offset%)
177 dels sout)
178 1.0f0))
179 (go label130)
180 label140
181 (cond
182 ((> lcode 2)
183 (setf iflag 6)
184 (go end_label)))
185 (f2cl-lib:fdo (jw 1 (f2cl-lib:int-add jw 1))
186 ((> jw np1) nil)
187 (tagbody
188 (setf (f2cl-lib:fref work-%data%
189 ((f2cl-lib:int-sub
190 (f2cl-lib:int-add iw jw)
193 (f2cl-lib:int-add
194 (f2cl-lib:int-mul 13
195 (f2cl-lib:int-add n
197 (f2cl-lib:int-mul 2 n)
198 lenqr)))
199 work-%offset%)
200 (qofs
201 (f2cl-lib:fref yold-%data%
202 (jw)
203 ((1 (f2cl-lib:int-add n 1)))
204 yold-%offset%)
205 (f2cl-lib:fref ypold-%data%
206 (jw)
207 ((1 (f2cl-lib:int-add n 1)))
208 ypold-%offset%)
209 (f2cl-lib:fref y-%data%
210 (jw)
211 ((1 (f2cl-lib:int-add n 1)))
212 y-%offset%)
213 (f2cl-lib:fref yp-%data%
214 (jw)
215 ((1 (f2cl-lib:int-add n 1)))
216 yp-%offset%)
217 dels sa))
218 label150))
219 (multiple-value-bind
220 (var-0 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-15 var-16)
222 (tangns sa
223 (f2cl-lib:array-slice work-%data%
224 double-float
225 (iw)
227 (f2cl-lib:int-add
228 (f2cl-lib:int-mul 13
229 (f2cl-lib:int-add
232 (f2cl-lib:int-mul 2 n)
233 lenqr)))
234 work-%offset%)
235 (f2cl-lib:array-slice work-%data%
236 double-float
237 (iwp)
239 (f2cl-lib:int-add
240 (f2cl-lib:int-mul 13
241 (f2cl-lib:int-add
244 (f2cl-lib:int-mul 2 n)
245 lenqr)))
246 work-%offset%)
247 (f2cl-lib:array-slice work-%data%
248 double-float
249 (itz)
251 (f2cl-lib:int-add
252 (f2cl-lib:int-mul 13
253 (f2cl-lib:int-add
256 (f2cl-lib:int-mul 2 n)
257 lenqr)))
258 work-%offset%)
259 ypold a qr lenqr pivot
260 (f2cl-lib:array-slice work-%data%
261 double-float
262 (ipp)
264 (f2cl-lib:int-add
265 (f2cl-lib:int-mul 13
266 (f2cl-lib:int-add
269 (f2cl-lib:int-mul 2 n)
270 lenqr)))
271 work-%offset%)
272 (f2cl-lib:array-slice work-%data%
273 double-float
274 (irho)
276 (f2cl-lib:int-add
277 (f2cl-lib:int-mul 13
278 (f2cl-lib:int-add
281 (f2cl-lib:int-mul 2 n)
282 lenqr)))
283 work-%offset%)
284 (f2cl-lib:array-slice work-%data%
285 double-float
286 (itangw)
288 (f2cl-lib:int-add
289 (f2cl-lib:int-mul 13
290 (f2cl-lib:int-add
293 (f2cl-lib:int-mul 2 n)
294 lenqr)))
295 work-%offset%)
296 nfe n iflag par ipar)
297 (declare (ignore var-1 var-2 var-3 var-4 var-5 var-6 var-7
298 var-8 var-9 var-10 var-11 var-13 var-15
299 var-16))
300 (setf sa var-0)
301 (setf nfe var-12)
302 (setf iflag var-14))
303 (if (> iflag 0) (go end_label))
304 (daxpy np1 1.0
305 (f2cl-lib:array-slice work-%data%
306 double-float
307 (itz)
309 (f2cl-lib:int-add
310 (f2cl-lib:int-mul 13
311 (f2cl-lib:int-add n
313 (f2cl-lib:int-mul 2 n)
314 lenqr)))
315 work-%offset%)
317 (f2cl-lib:array-slice work-%data%
318 double-float
319 (iw)
321 (f2cl-lib:int-add
322 (f2cl-lib:int-mul 13
323 (f2cl-lib:int-add n
325 (f2cl-lib:int-mul 2 n)
326 lenqr)))
327 work-%offset%)
329 (multiple-value-bind
330 (var-0 var-1 var-2 var-3 var-4 var-5 var-6 var-7 var-8
331 var-9 var-10 var-11 var-12 var-13 var-14 var-15 var-16)
332 (tangns sa
333 (f2cl-lib:array-slice work-%data%
334 double-float
335 (iw)
337 (f2cl-lib:int-add
338 (f2cl-lib:int-mul 13
339 (f2cl-lib:int-add
342 (f2cl-lib:int-mul 2 n)
343 lenqr)))
344 work-%offset%)
345 (f2cl-lib:array-slice work-%data%
346 double-float
347 (iwp)
349 (f2cl-lib:int-add
350 (f2cl-lib:int-mul 13
351 (f2cl-lib:int-add
354 (f2cl-lib:int-mul 2 n)
355 lenqr)))
356 work-%offset%)
357 (f2cl-lib:array-slice work-%data%
358 double-float
359 (itz)
361 (f2cl-lib:int-add
362 (f2cl-lib:int-mul 13
363 (f2cl-lib:int-add
366 (f2cl-lib:int-mul 2 n)
367 lenqr)))
368 work-%offset%)
369 ypold a qr lenqr pivot
370 (f2cl-lib:array-slice work-%data%
371 double-float
372 (ipp)
374 (f2cl-lib:int-add
375 (f2cl-lib:int-mul 13
376 (f2cl-lib:int-add
379 (f2cl-lib:int-mul 2 n)
380 lenqr)))
381 work-%offset%)
382 (f2cl-lib:array-slice work-%data%
383 double-float
384 (irho)
386 (f2cl-lib:int-add
387 (f2cl-lib:int-mul 13
388 (f2cl-lib:int-add
391 (f2cl-lib:int-mul 2 n)
392 lenqr)))
393 work-%offset%)
394 (f2cl-lib:array-slice work-%data%
395 double-float
396 (itangw)
398 (f2cl-lib:int-add
399 (f2cl-lib:int-mul 13
400 (f2cl-lib:int-add
403 (f2cl-lib:int-mul 2 n)
404 lenqr)))
405 work-%offset%)
406 nfe n iflag par ipar)
407 (declare (ignore var-1 var-2 var-3 var-4 var-5 var-6 var-7
408 var-8 var-9 var-10 var-11 var-13 var-15
409 var-16))
410 (setf sa var-0)
411 (setf nfe var-12)
412 (setf iflag var-14))
413 (if (> iflag 0) (go end_label))
414 (daxpy np1 1.0
415 (f2cl-lib:array-slice work-%data%
416 double-float
417 (itz)
419 (f2cl-lib:int-add
420 (f2cl-lib:int-mul 13
421 (f2cl-lib:int-add n
423 (f2cl-lib:int-mul 2 n)
424 lenqr)))
425 work-%offset%)
427 (f2cl-lib:array-slice work-%data%
428 double-float
429 (iw)
431 (f2cl-lib:int-add
432 (f2cl-lib:int-mul 13
433 (f2cl-lib:int-add n
435 (f2cl-lib:int-mul 2 n)
436 lenqr)))
437 work-%offset%)
439 (cond
440 ((and
442 (abs
444 (f2cl-lib:fref work
445 ((f2cl-lib:int-add iw n))
447 (f2cl-lib:int-add
448 (f2cl-lib:int-mul 13
449 (f2cl-lib:int-add n
451 (f2cl-lib:int-mul 2 n)
452 lenqr))))
453 (- 1.0f0)))
454 (+ rerr aerr))
456 (dnrm2 n
457 (f2cl-lib:array-slice work
458 double-float
459 (itz)
461 (f2cl-lib:int-add
462 (f2cl-lib:int-mul 13
463 (f2cl-lib:int-add
466 (f2cl-lib:int-mul 2 n)
467 lenqr))))
470 (* rerr
471 (dnrm2 n
472 (f2cl-lib:array-slice work
473 double-float
474 (iw)
476 (f2cl-lib:int-add
477 (f2cl-lib:int-mul 13
478 (f2cl-lib:int-add
481 (f2cl-lib:int-mul 2 n)
482 lenqr))))
484 aerr)))
485 (dcopy np1
486 (f2cl-lib:array-slice work-%data%
487 double-float
488 (iw)
490 (f2cl-lib:int-add
491 (f2cl-lib:int-mul 13
492 (f2cl-lib:int-add
495 (f2cl-lib:int-mul 2 n)
496 lenqr)))
497 work-%offset%)
498 1 y 1)
499 (go end_label)))
500 (cond
503 (+ (f2cl-lib:fref yold (np1) ((1 (f2cl-lib:int-add n 1))))
504 (- 1.0f0))
506 (f2cl-lib:fref work
507 ((f2cl-lib:int-add iw n))
509 (f2cl-lib:int-add
510 (f2cl-lib:int-mul 13
511 (f2cl-lib:int-add n 1))
512 (f2cl-lib:int-mul 2 n)
513 lenqr))))
514 (- 1.0f0)))
515 0.0f0)
516 (dcopy np1
517 (f2cl-lib:array-slice work-%data%
518 double-float
519 (iw)
521 (f2cl-lib:int-add
522 (f2cl-lib:int-mul 13
523 (f2cl-lib:int-add
526 (f2cl-lib:int-mul 2 n)
527 lenqr)))
528 work-%offset%)
529 1 yold 1)
530 (dcopy np1
531 (f2cl-lib:array-slice work-%data%
532 double-float
533 (iwp)
535 (f2cl-lib:int-add
536 (f2cl-lib:int-mul 13
537 (f2cl-lib:int-add
540 (f2cl-lib:int-mul 2 n)
541 lenqr)))
542 work-%offset%)
543 1 ypold 1))
545 (dcopy np1
546 (f2cl-lib:array-slice work-%data%
547 double-float
548 (iw)
550 (f2cl-lib:int-add
551 (f2cl-lib:int-mul 13
552 (f2cl-lib:int-add
555 (f2cl-lib:int-mul 2 n)
556 lenqr)))
557 work-%offset%)
558 1 y 1)
559 (dcopy np1
560 (f2cl-lib:array-slice work-%data%
561 double-float
562 (iwp)
564 (f2cl-lib:int-add
565 (f2cl-lib:int-mul 13
566 (f2cl-lib:int-add
569 (f2cl-lib:int-mul 2 n)
570 lenqr)))
571 work-%offset%)
572 1 yp 1)))
573 label300))
574 (setf iflag 6)
575 (go end_label)
576 end_label
577 (return
578 (values nil
580 iflag
593 nil)))))))
595 (in-package #-gcl #:cl-user #+gcl "CL-USER")
596 #+#.(cl:if (cl:find-package '#:f2cl) '(and) '(or))
597 (eval-when (:load-toplevel :compile-toplevel :execute)
598 (setf (gethash 'fortran-to-lisp::rootns
599 fortran-to-lisp::*f2cl-function-info*)
600 (fortran-to-lisp::make-f2cl-finfo
601 :arg-types '((fortran-to-lisp::integer4) (fortran-to-lisp::integer4)
602 (fortran-to-lisp::integer4) (double-float)
603 (double-float) (array double-float (*))
604 (array double-float (*)) (array double-float (*))
605 (array double-float (*)) (array double-float (*))
606 (array double-float (*)) (fortran-to-lisp::integer4)
607 (array fortran-to-lisp::integer4 (*))
608 (array double-float (*)) (array double-float (*))
609 (array fortran-to-lisp::integer4 (*)))
610 :return-values '(nil fortran-to-lisp::nfe fortran-to-lisp::iflag nil
611 nil nil nil nil nil nil nil nil nil nil nil nil)
612 :calls '(fortran-to-lisp::dcopy fortran-to-lisp::daxpy
613 fortran-to-lisp::dnrm2 fortran-to-lisp::tangns
614 fortran-to-lisp::root fortran-to-lisp::d1mach))))