Remove the obsolete DEFMTRFUN-EXTERNAL macro
[maxima.git] / share / hompack / lisp / tangnf.lisp
blob36ed5f4cf469dbe8bebdab63868b142c85cc7446
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 tangnf (rholen y yp ypold a qr alpha tz pivot nfe n iflag par ipar)
21 (declare (type (array f2cl-lib:integer4 (*)) ipar)
22 (type (array double-float (*)) par)
23 (type (f2cl-lib:integer4) iflag n nfe)
24 (type (array f2cl-lib:integer4 (*)) pivot)
25 (type (array double-float (*)) tz alpha qr a ypold yp y)
26 (type (double-float) rholen))
27 (f2cl-lib:with-multi-array-data
28 ((y double-float y-%data% y-%offset%)
29 (yp double-float yp-%data% yp-%offset%)
30 (ypold double-float ypold-%data% ypold-%offset%)
31 (a double-float a-%data% a-%offset%)
32 (qr double-float qr-%data% qr-%offset%)
33 (alpha double-float alpha-%data% alpha-%offset%)
34 (tz double-float tz-%data% tz-%offset%)
35 (pivot f2cl-lib:integer4 pivot-%data% pivot-%offset%)
36 (par double-float par-%data% par-%offset%)
37 (ipar f2cl-lib:integer4 ipar-%data% ipar-%offset%))
38 (prog ((i 0) (j 0) (jbar 0) (k 0) (kp1 0) (np1 0) (np2 0) (alphak 0.0)
39 (beta 0.0) (qrkk 0.0) (sigma 0.0) (sum 0.0) (ypnorm 0.0)
40 (lambda$ 0.0))
41 (declare (type (double-float) lambda$ ypnorm sum sigma qrkk beta alphak)
42 (type (f2cl-lib:integer4) np2 np1 kp1 k jbar j i))
43 (setf lambda$
44 (f2cl-lib:fref y-%data%
45 (1)
46 ((1 (f2cl-lib:int-add n 1)))
47 y-%offset%))
48 (setf np1 (f2cl-lib:int-add n 1))
49 (setf np2 (f2cl-lib:int-add n 2))
50 (setf nfe (f2cl-lib:int-add nfe 1))
51 (cond
52 ((= iflag (f2cl-lib:int-sub 2))
53 (f2cl-lib:fdo (k 1 (f2cl-lib:int-add k 1))
54 ((> k np1) nil)
55 (tagbody
56 (multiple-value-bind (var-0 var-1 var-2 var-3 var-4 var-5 var-6)
57 (rhojac a lambda$
58 (f2cl-lib:array-slice y-%data%
59 double-float
60 (2)
61 ((1 (f2cl-lib:int-add n 1)))
62 y-%offset%)
63 (f2cl-lib:array-slice qr-%data%
64 double-float
65 (1 k)
66 ((1 n) (1 (f2cl-lib:int-add n 2)))
67 qr-%offset%)
68 k par ipar)
69 (declare (ignore var-0 var-2 var-3 var-4 var-5 var-6))
70 (setf lambda$ var-1))
71 label30))
72 (multiple-value-bind (var-0 var-1 var-2 var-3 var-4 var-5)
73 (rho a lambda$
74 (f2cl-lib:array-slice y-%data%
75 double-float
76 (2)
77 ((1 (f2cl-lib:int-add n 1)))
78 y-%offset%)
79 (f2cl-lib:array-slice qr-%data%
80 double-float
81 (1 np2)
82 ((1 n) (1 (f2cl-lib:int-add n 2)))
83 qr-%offset%)
84 par ipar)
85 (declare (ignore var-0 var-2 var-3 var-4 var-5))
86 (setf lambda$ var-1)))
89 (f2cl-lib:array-slice y-%data%
90 double-float
91 (2)
92 ((1 (f2cl-lib:int-add n 1)))
93 y-%offset%)
94 tz)
95 (cond
96 ((= iflag 0)
97 (f2cl-lib:fdo (j 1 (f2cl-lib:int-add j 1))
98 ((> j n) nil)
99 (tagbody
100 (setf sigma (f2cl-lib:fref a-%data% (j) ((1 n)) a-%offset%))
101 (setf beta
102 (- sigma
103 (f2cl-lib:fref tz-%data%
105 ((1 (f2cl-lib:int-add n 1)))
106 tz-%offset%)))
107 (setf (f2cl-lib:fref qr-%data%
108 (j 1)
109 ((1 n) (1 (f2cl-lib:int-add n 2)))
110 qr-%offset%)
111 beta)
112 label100
113 (setf (f2cl-lib:fref qr-%data%
114 (j np2)
115 ((1 n) (1 (f2cl-lib:int-add n 2)))
116 qr-%offset%)
119 (f2cl-lib:fref y-%data%
120 ((f2cl-lib:int-add j 1))
121 ((1 (f2cl-lib:int-add n 1)))
122 y-%offset%)
123 sigma)
124 (* lambda$ beta)))))
125 (f2cl-lib:fdo (k 1 (f2cl-lib:int-add k 1))
126 ((> k n) nil)
127 (tagbody
128 (fjac
129 (f2cl-lib:array-slice y-%data%
130 double-float
132 ((1 (f2cl-lib:int-add n 1)))
133 y-%offset%)
134 tz k)
135 (setf kp1 (f2cl-lib:int-add k 1))
136 (f2cl-lib:fdo (j 1 (f2cl-lib:int-add j 1))
137 ((> j n) nil)
138 (tagbody
139 label110
140 (setf (f2cl-lib:fref qr-%data%
141 (j kp1)
142 ((1 n) (1 (f2cl-lib:int-add n 2)))
143 qr-%offset%)
144 (* (- lambda$)
145 (f2cl-lib:fref tz-%data%
147 ((1 (f2cl-lib:int-add n 1)))
148 tz-%offset%)))))
149 label120
150 (setf (f2cl-lib:fref qr-%data%
151 (k kp1)
152 ((1 n) (1 (f2cl-lib:int-add n 2)))
153 qr-%offset%)
154 (+ 1.0f0
155 (f2cl-lib:fref qr-%data%
156 (k kp1)
157 ((1 n) (1 (f2cl-lib:int-add n 2)))
158 qr-%offset%))))))
160 (tagbody
161 label140
162 (f2cl-lib:fdo (j 1 (f2cl-lib:int-add j 1))
163 ((> j n) nil)
164 (tagbody
165 (setf sigma
167 (f2cl-lib:fref y-%data%
168 ((f2cl-lib:int-add j 1))
169 ((1 (f2cl-lib:int-add n 1)))
170 y-%offset%)
171 (f2cl-lib:fref a-%data% (j) ((1 n)) a-%offset%)))
172 (setf beta
174 (f2cl-lib:fref tz-%data%
176 ((1 (f2cl-lib:int-add n 1)))
177 tz-%offset%)
178 sigma))
179 (setf (f2cl-lib:fref qr-%data%
180 (j 1)
181 ((1 n) (1 (f2cl-lib:int-add n 2)))
182 qr-%offset%)
183 beta)
184 label150
185 (setf (f2cl-lib:fref qr-%data%
186 (j np2)
187 ((1 n) (1 (f2cl-lib:int-add n 2)))
188 qr-%offset%)
189 (+ sigma (* lambda$ beta)))))
190 (f2cl-lib:fdo (k 1 (f2cl-lib:int-add k 1))
191 ((> k n) nil)
192 (tagbody
193 (fjac
194 (f2cl-lib:array-slice y-%data%
195 double-float
197 ((1 (f2cl-lib:int-add n 1)))
198 y-%offset%)
199 tz k)
200 (setf kp1 (f2cl-lib:int-add k 1))
201 (f2cl-lib:fdo (j 1 (f2cl-lib:int-add j 1))
202 ((> j n) nil)
203 (tagbody
204 label160
205 (setf (f2cl-lib:fref qr-%data%
206 (j kp1)
207 ((1 n) (1 (f2cl-lib:int-add n 2)))
208 qr-%offset%)
209 (* lambda$
210 (f2cl-lib:fref tz-%data%
212 ((1 (f2cl-lib:int-add n 1)))
213 tz-%offset%)))))
214 label170
215 (setf (f2cl-lib:fref qr-%data%
216 (k kp1)
217 ((1 n) (1 (f2cl-lib:int-add n 2)))
218 qr-%offset%)
219 (+ (- 1.0f0 lambda$)
220 (f2cl-lib:fref qr-%data%
221 (k kp1)
222 ((1 n) (1 (f2cl-lib:int-add n 2)))
223 qr-%offset%))))))))))
224 (if (< rholen 0.0f0)
225 (setf rholen
226 (dnrm2 n
227 (f2cl-lib:array-slice qr-%data%
228 double-float
229 (1 np2)
230 ((1 n) (1 (f2cl-lib:int-add n 2)))
231 qr-%offset%)
232 1)))
233 (f2cl-lib:fdo (j 1 (f2cl-lib:int-add j 1))
234 ((> j np1) nil)
235 (tagbody
236 (setf (f2cl-lib:fref yp-%data%
238 ((1 (f2cl-lib:int-add n 1)))
239 yp-%offset%)
240 (ddot n
241 (f2cl-lib:array-slice qr-%data%
242 double-float
243 (1 j)
244 ((1 n) (1 (f2cl-lib:int-add n 2)))
245 qr-%offset%)
247 (f2cl-lib:array-slice qr-%data%
248 double-float
249 (1 j)
250 ((1 n) (1 (f2cl-lib:int-add n 2)))
251 qr-%offset%)
253 label220
254 (setf (f2cl-lib:fref pivot-%data%
256 ((1 (f2cl-lib:int-add n 1)))
257 pivot-%offset%)
258 j)))
259 (f2cl-lib:fdo (k 1 (f2cl-lib:int-add k 1))
260 ((> k n) nil)
261 (tagbody
262 (setf sigma
263 (f2cl-lib:fref yp-%data%
265 ((1 (f2cl-lib:int-add n 1)))
266 yp-%offset%))
267 (setf jbar k)
268 (setf kp1 (f2cl-lib:int-add k 1))
269 (f2cl-lib:fdo (j kp1 (f2cl-lib:int-add j 1))
270 ((> j np1) nil)
271 (tagbody
273 (>= sigma
274 (f2cl-lib:fref yp-%data%
276 ((1 (f2cl-lib:int-add n 1)))
277 yp-%offset%))
278 (go label240))
279 (setf sigma
280 (f2cl-lib:fref yp-%data%
282 ((1 (f2cl-lib:int-add n 1)))
283 yp-%offset%))
284 (setf jbar j)
285 label240))
286 (if (= jbar k) (go label260))
287 (setf i
288 (f2cl-lib:fref pivot-%data%
290 ((1 (f2cl-lib:int-add n 1)))
291 pivot-%offset%))
292 (setf (f2cl-lib:fref pivot-%data%
294 ((1 (f2cl-lib:int-add n 1)))
295 pivot-%offset%)
296 (f2cl-lib:fref pivot-%data%
297 (jbar)
298 ((1 (f2cl-lib:int-add n 1)))
299 pivot-%offset%))
300 (setf (f2cl-lib:fref pivot-%data%
301 (jbar)
302 ((1 (f2cl-lib:int-add n 1)))
303 pivot-%offset%)
305 (setf (f2cl-lib:fref yp-%data%
306 (jbar)
307 ((1 (f2cl-lib:int-add n 1)))
308 yp-%offset%)
309 (f2cl-lib:fref yp-%data%
311 ((1 (f2cl-lib:int-add n 1)))
312 yp-%offset%))
313 (setf (f2cl-lib:fref yp-%data%
315 ((1 (f2cl-lib:int-add n 1)))
316 yp-%offset%)
317 sigma)
318 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1))
319 ((> i n) nil)
320 (tagbody
321 (setf sigma
322 (f2cl-lib:fref qr-%data%
323 (i k)
324 ((1 n) (1 (f2cl-lib:int-add n 2)))
325 qr-%offset%))
326 (setf (f2cl-lib:fref qr-%data%
327 (i k)
328 ((1 n) (1 (f2cl-lib:int-add n 2)))
329 qr-%offset%)
330 (f2cl-lib:fref qr-%data%
331 (i jbar)
332 ((1 n) (1 (f2cl-lib:int-add n 2)))
333 qr-%offset%))
334 (setf (f2cl-lib:fref qr-%data%
335 (i jbar)
336 ((1 n) (1 (f2cl-lib:int-add n 2)))
337 qr-%offset%)
338 sigma)
339 label250))
340 label260
341 (setf sigma
342 (ddot (f2cl-lib:int-add (f2cl-lib:int-sub n k) 1)
343 (f2cl-lib:array-slice qr-%data%
344 double-float
345 (k k)
346 ((1 n) (1 (f2cl-lib:int-add n 2)))
347 qr-%offset%)
349 (f2cl-lib:array-slice qr-%data%
350 double-float
351 (k k)
352 ((1 n) (1 (f2cl-lib:int-add n 2)))
353 qr-%offset%)
355 (cond
356 ((= sigma 0.0f0)
357 (setf iflag 4)
358 (go end_label)))
359 label270
360 (if (= k n) (go label300))
361 (setf qrkk
362 (f2cl-lib:fref qr-%data%
363 (k k)
364 ((1 n) (1 (f2cl-lib:int-add n 2)))
365 qr-%offset%))
366 (setf alphak (- (f2cl-lib:fsqrt sigma)))
367 (if (< qrkk 0.0f0) (setf alphak (- alphak)))
368 (setf (f2cl-lib:fref alpha-%data% (k) ((1 n)) alpha-%offset%) alphak)
369 (setf beta (/ 1.0f0 (- sigma (* qrkk alphak))))
370 (setf (f2cl-lib:fref qr-%data%
371 (k k)
372 ((1 n) (1 (f2cl-lib:int-add n 2)))
373 qr-%offset%)
374 (- qrkk alphak))
375 (f2cl-lib:fdo (j kp1 (f2cl-lib:int-add j 1))
376 ((> j np2) nil)
377 (tagbody
378 (setf sigma
379 (* beta
380 (ddot (f2cl-lib:int-add (f2cl-lib:int-sub n k) 1)
381 (f2cl-lib:array-slice qr-%data%
382 double-float
383 (k k)
384 ((1 n)
385 (1 (f2cl-lib:int-add n 2)))
386 qr-%offset%)
388 (f2cl-lib:array-slice qr-%data%
389 double-float
390 (k j)
391 ((1 n)
392 (1 (f2cl-lib:int-add n 2)))
393 qr-%offset%)
394 1)))
395 (f2cl-lib:fdo (i k (f2cl-lib:int-add i 1))
396 ((> i n) nil)
397 (tagbody
398 (setf (f2cl-lib:fref qr-%data%
399 (i j)
400 ((1 n) (1 (f2cl-lib:int-add n 2)))
401 qr-%offset%)
403 (f2cl-lib:fref qr-%data%
404 (i j)
405 ((1 n) (1 (f2cl-lib:int-add n 2)))
406 qr-%offset%)
408 (f2cl-lib:fref qr-%data%
409 (i k)
410 ((1 n) (1 (f2cl-lib:int-add n 2)))
411 qr-%offset%)
412 sigma)))
413 label280))
414 (if (< j np2)
415 (setf (f2cl-lib:fref yp-%data%
417 ((1 (f2cl-lib:int-add n 1)))
418 yp-%offset%)
420 (f2cl-lib:fref yp-%data%
422 ((1 (f2cl-lib:int-add n 1)))
423 yp-%offset%)
424 (expt
425 (f2cl-lib:fref qr-%data%
426 (k j)
427 ((1 n) (1 (f2cl-lib:int-add n 2)))
428 qr-%offset%)
429 2))))
430 label290))
431 label300))
432 (setf (f2cl-lib:fref alpha-%data% (n) ((1 n)) alpha-%offset%)
433 (f2cl-lib:fref qr-%data%
434 (n n)
435 ((1 n) (1 (f2cl-lib:int-add n 2)))
436 qr-%offset%))
437 (setf (f2cl-lib:fref tz-%data%
438 (np1)
439 ((1 (f2cl-lib:int-add n 1)))
440 tz-%offset%)
441 (coerce 1.0f0 'double-float))
442 (f2cl-lib:fdo (i n (f2cl-lib:int-add i (f2cl-lib:int-sub 1)))
443 ((> i 1) nil)
444 (tagbody
445 (setf sum (coerce 0.0f0 'double-float))
446 (f2cl-lib:fdo (j (f2cl-lib:int-add i 1) (f2cl-lib:int-add j 1))
447 ((> j np1) nil)
448 (tagbody
449 label330
450 (setf sum
451 (+ sum
453 (f2cl-lib:fref qr-%data%
454 (i j)
455 ((1 n) (1 (f2cl-lib:int-add n 2)))
456 qr-%offset%)
457 (f2cl-lib:fref tz-%data%
459 ((1 (f2cl-lib:int-add n 1)))
460 tz-%offset%))))))
461 label340
462 (setf (f2cl-lib:fref tz-%data%
464 ((1 (f2cl-lib:int-add n 1)))
465 tz-%offset%)
466 (/ (- sum)
467 (f2cl-lib:fref alpha-%data%
469 ((1 n))
470 alpha-%offset%)))))
471 (setf ypnorm (dnrm2 np1 tz 1))
472 (f2cl-lib:fdo (k 1 (f2cl-lib:int-add k 1))
473 ((> k np1) nil)
474 (tagbody
475 label360
476 (setf (f2cl-lib:fref yp-%data%
477 ((f2cl-lib:fref pivot
479 ((1 (f2cl-lib:int-add n 1)))))
480 ((1 (f2cl-lib:int-add n 1)))
481 yp-%offset%)
483 (f2cl-lib:fref tz-%data%
485 ((1 (f2cl-lib:int-add n 1)))
486 tz-%offset%)
487 ypnorm))))
488 (if (>= (ddot np1 yp 1 ypold 1) 0.0f0) (go label380))
489 (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1))
490 ((> i np1) nil)
491 (tagbody
492 label370
493 (setf (f2cl-lib:fref yp-%data%
495 ((1 (f2cl-lib:int-add n 1)))
496 yp-%offset%)
498 (f2cl-lib:fref yp-%data%
500 ((1 (f2cl-lib:int-add n 1)))
501 yp-%offset%)))))
502 label380
503 (f2cl-lib:fdo (i n (f2cl-lib:int-add i (f2cl-lib:int-sub 1)))
504 ((> i 1) nil)
505 (tagbody
506 (setf sum
508 (f2cl-lib:fref qr-%data%
509 (i np1)
510 ((1 n) (1 (f2cl-lib:int-add n 2)))
511 qr-%offset%)
512 (f2cl-lib:fref qr-%data%
513 (i np2)
514 ((1 n) (1 (f2cl-lib:int-add n 2)))
515 qr-%offset%)))
516 (f2cl-lib:fdo (j (f2cl-lib:int-add i 1) (f2cl-lib:int-add j 1))
517 ((> j n) nil)
518 (tagbody
519 label430
520 (setf sum
521 (+ sum
523 (f2cl-lib:fref qr-%data%
524 (i j)
525 ((1 n) (1 (f2cl-lib:int-add n 2)))
526 qr-%offset%)
527 (f2cl-lib:fref alpha-%data%
529 ((1 n))
530 alpha-%offset%))))))
531 label440
532 (setf (f2cl-lib:fref alpha-%data% (i) ((1 n)) alpha-%offset%)
533 (/ (- sum)
534 (f2cl-lib:fref alpha-%data%
536 ((1 n))
537 alpha-%offset%)))))
538 (f2cl-lib:fdo (k 1 (f2cl-lib:int-add k 1))
539 ((> k n) nil)
540 (tagbody
541 label450
542 (setf (f2cl-lib:fref tz-%data%
543 ((f2cl-lib:fref pivot
545 ((1 (f2cl-lib:int-add n 1)))))
546 ((1 (f2cl-lib:int-add n 1)))
547 tz-%offset%)
548 (f2cl-lib:fref alpha-%data% (k) ((1 n)) alpha-%offset%))))
549 (setf (f2cl-lib:fref tz-%data%
550 ((f2cl-lib:fref pivot
551 (np1)
552 ((1 (f2cl-lib:int-add n 1)))))
553 ((1 (f2cl-lib:int-add n 1)))
554 tz-%offset%)
555 (coerce 1.0f0 'double-float))
556 (setf sigma (ddot np1 tz 1 yp 1))
557 (f2cl-lib:fdo (j 1 (f2cl-lib:int-add j 1))
558 ((> j np1) nil)
559 (tagbody
560 (setf (f2cl-lib:fref tz-%data%
562 ((1 (f2cl-lib:int-add n 1)))
563 tz-%offset%)
565 (f2cl-lib:fref tz-%data%
567 ((1 (f2cl-lib:int-add n 1)))
568 tz-%offset%)
569 (* sigma
570 (f2cl-lib:fref yp-%data%
572 ((1 (f2cl-lib:int-add n 1)))
573 yp-%offset%))))
574 label470))
575 (go end_label)
576 end_label
577 (return
578 (values rholen nil nil nil nil nil nil nil nil nfe nil iflag nil nil)))))
580 (in-package #:cl-user)
581 #+#.(cl:if (cl:find-package '#:f2cl) '(and) '(or))
582 (eval-when (:load-toplevel :compile-toplevel :execute)
583 (setf (gethash 'fortran-to-lisp::tangnf
584 fortran-to-lisp::*f2cl-function-info*)
585 (fortran-to-lisp::make-f2cl-finfo
586 :arg-types '((double-float) (array double-float (*))
587 (array double-float (*)) (array double-float (*))
588 (array double-float (*)) (array double-float (*))
589 (array double-float (*)) (array double-float (*))
590 (array fortran-to-lisp::integer4 (*))
591 (fortran-to-lisp::integer4) (fortran-to-lisp::integer4)
592 (fortran-to-lisp::integer4) (array double-float (*))
593 (array fortran-to-lisp::integer4 (*)))
594 :return-values '(fortran-to-lisp::rholen nil nil nil nil nil nil nil
595 nil fortran-to-lisp::nfe nil fortran-to-lisp::iflag
596 nil nil)
597 :calls '(fortran-to-lisp::fjac fortran-to-lisp::f
598 fortran-to-lisp::ddot fortran-to-lisp::dnrm2
599 fortran-to-lisp::rho fortran-to-lisp::rhojac))))