linalg 1.6.1
A linear algebra library that provides a user-friendly interface to several BLAS and LAPACK routines.
Loading...
Searching...
No Matches
linalg_c_api.f90
1! linalg_c_api.f90
2
6 use iso_c_binding
7 use linalg
8 use ferror
9 implicit none
10
11contains
12! ------------------------------------------------------------------------------
29 function la_rank1_update(m, n, alpha, x, y, a, lda) &
30 bind(C, name = "la_rank1_update") result(flag)
31 ! Arguments
32 integer(c_int), intent(in), value :: m, n, lda
33 real(c_double), intent(in), value :: alpha
34 real(c_double), intent(in) :: x(*), y(*)
35 real(c_double), intent(inout) :: a(lda,*)
36 integer(c_int) :: flag
37
38 ! Initialization
39 flag = la_no_error
40
41 ! Input Checking
42 if (lda < m) then
43 flag = la_invalid_input_error
44 return
45 end if
46
47 ! Process
48 call rank1_update(alpha, x(1:m), y(1:n), a(1:m,1:n))
49 end function
50
51! ------------------------------------------------------------------------------
68 function la_rank1_update_cmplx(m, n, alpha, x, y, a, lda) &
69 bind(C, name = "la_rank1_update_cmplx") result(flag)
70 ! Arguments
71 integer(c_int), intent(in), value :: m, n, lda
72 complex(c_double), intent(in), value :: alpha
73 complex(c_double), intent(in) :: x(*), y(*)
74 complex(c_double), intent(inout) :: a(lda,*)
75 integer(c_int) :: flag
76
77 ! Initialization
78 flag = la_no_error
79
80 ! Input Checking
81 if (lda < m) then
82 flag = la_invalid_input_error
83 return
84 end if
85
86 ! Process
87 call rank1_update(alpha, x(1:m), y(1:n), a(1:m,1:n))
88 end function
89
90! ------------------------------------------------------------------------------
103 function la_trace(m, n, a, lda, rst) bind(C, name = "la_trace") &
104 result(flag)
105 ! Arguments
106 integer(c_int), intent(in), value :: m, n, lda
107 real(c_double), intent(in) :: a(lda,*)
108 real(c_double), intent(out) :: rst
109 integer(c_int) :: flag
110
111 ! Initialization
112 flag = la_no_error
113
114 ! Input Checking
115 if (lda < m) then
116 flag = la_invalid_input_error
117 return
118 end if
119
120 ! Process
121 rst = trace(a(1:m,1:n))
122 end function
123
124! ------------------------------------------------------------------------------
137 function la_trace_cmplx(m, n, a, lda, rst) &
138 bind(C, name = "la_trace_cmplx") result(flag)
139 ! Arguments
140 integer(c_int), intent(in), value :: m, n, lda
141 complex(c_double), intent(in) :: a(lda,*)
142 complex(c_double), intent(out) :: rst
143 integer(c_int) :: flag
144
145 ! Initialization
146 flag = la_no_error
147
148 ! Input Checking
149 if (lda < m) then
150 flag = la_invalid_input_error
151 return
152 end if
153
154 ! Process
155 rst = trace(a(1:m,1:n))
156 end function
157
158! ------------------------------------------------------------------------------
183 function la_mtx_mult(transa, transb, m, n, k, alpha, a, lda, b, ldb, &
184 beta, c, ldc) bind(C, name="la_mtx_mult") result(flag)
185 ! Arugments
186 logical(c_bool), intent(in), value :: transa, transb
187 integer(c_int), intent(in), value :: m, n, k, lda, ldb, ldc
188 real(c_double), intent(in), value :: alpha, beta
189 real(c_double), intent(in) :: a(lda,*), b(ldb,*)
190 real(c_double), intent(inout) :: c(ldc,*)
191 integer(c_int) :: flag
192
193 ! Local Variables
194 character :: ta, tb
195 integer(c_int) :: nrowa, nrowb
196
197 ! Initialization
198 flag = la_no_error
199 ta = "N"
200 if (transa) ta = "T"
201
202 tb = "N"
203 if (transb) tb = "T"
204
205 if (transa) then
206 nrowa = k
207 else
208 nrowa = m
209 end if
210
211 if (transb) then
212 nrowb = n
213 else
214 nrowb = k
215 end if
216
217 ! Input Checking
218 if (lda < nrowa .or. ldb < nrowb .or. ldc < m) then
219 flag = la_invalid_input_error
220 return
221 end if
222
223 ! Call DGEMM directly
224 call dgemm(ta, tb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
225 end function
226
227! ------------------------------------------------------------------------------
254 function la_mtx_mult_cmplx(opa, opb, m, n, k, alpha, a, lda, b, ldb, &
255 beta, c, ldc) bind(C, name="la_mtx_mult_cmplx") result(flag)
256 ! Arguments
257 integer(c_int), intent(in), value :: opa, opb, m, n, k, lda, ldb, ldc
258 complex(c_double), intent(in), value :: alpha, beta
259 complex(c_double), intent(in) :: a(lda,*), b(ldb,*)
260 complex(c_double), intent(inout) :: c(ldc,*)
261 integer(c_int) :: flag
262
263 ! Local Variables
264 character :: ta, tb
265 integer(c_int) :: nrowa, nrowb
266
267 ! Initialization
268 flag = la_no_error
269 if (opa == la_transpose) then
270 ta = "T"
271 else if (opa == la_hermitian_transpose) then
272 ta = "H"
273 else
274 ta = "N"
275 end if
276
277 if (opb == la_transpose) then
278 tb = "T"
279 else if (opb == la_hermitian_transpose) then
280 tb = "H"
281 else
282 tb = "N"
283 end if
284
285 if (opa == la_transpose .or. opa == la_hermitian_transpose) then
286 nrowa = k
287 else
288 nrowa = m
289 end if
290
291 if (opb == la_transpose .or. opb == la_hermitian_transpose) then
292 nrowb = n
293 else
294 nrowb = k
295 end if
296
297 ! Input Checking
298 if (lda < nrowa .or. ldb < nrowb .or. ldc < m) then
299 flag = la_invalid_input_error
300 return
301 end if
302
303 ! Call ZGEMM directly
304 call zgemm(ta, tb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
305 end function
306
307! ------------------------------------------------------------------------------
339 function la_diag_mtx_mult(lside, transb, m, n, k, alpha, a, b, ldb, &
340 beta, c, ldc) bind(C, name="la_diag_mtx_mult") result(flag)
341 ! Arguments
342 logical(c_bool), intent(in), value :: lside, transb
343 integer(c_int), intent(in), value :: m, n, k, ldb, ldc
344 real(c_double), intent(in), value :: alpha, beta
345 real(c_double), intent(in) :: a(*), b(ldb,*)
346 real(c_double), intent(inout) :: c(ldc,*)
347 integer(c_int) :: flag
348
349 ! Local Variabes
350 integer(c_int) :: nrows, ncols, p
351 logical :: ls, tb
352 type(errors) :: err
353
354 ! Initialization
355 call err%set_exit_on_error(.false.)
356 flag = la_no_error
357 if (lside .and. transb) then
358 nrows = n
359 ncols = k
360 p = min(k, m)
361 ls = .true.
362 tb = .true.
363 else if (lside .and. .not. transb) then
364 nrows = k
365 ncols = n
366 p = min(k, m)
367 ls = .true.
368 tb = .false.
369 else if (.not. lside .and. transb) then
370 nrows = k
371 ncols = m
372 p = min(k, n)
373 ls = .false.
374 tb = .true.
375 else
376 nrows = m
377 ncols = k
378 p = min(k, n)
379 ls = .false.
380 tb = .false.
381 end if
382
383 ! Error Checking
384 if (ldb < nrows .or. ldc < m) then
385 flag = la_invalid_input_error
386 return
387 end if
388
389 ! Process
390 call diag_mtx_mult(ls, tb, alpha, a(1:p), b(1:nrows,1:ncols), &
391 beta, c(1:m,1:n), err)
392 if (err%has_error_occurred()) flag = err%get_error_flag()
393 end function
394
395! ------------------------------------------------------------------------------
428 function la_diag_mtx_mult_mixed(lside, opb, m, n, k, alpha, a, b, ldb, &
429 beta, c, ldc) bind(C, name = "la_diag_mtx_mult_mixed") result(flag)
430 ! Arguments
431 logical(c_bool), intent(in), value :: lside
432 integer(c_int), intent(in), value :: opb, m, n, k, ldb, ldc
433 complex(c_double), intent(in), value :: alpha, beta
434 real(c_double), intent(in) :: a(*)
435 complex(c_double), intent(in) :: b(ldb,*)
436 complex(c_double), intent(inout) :: c(ldc,*)
437 integer(c_int) :: flag
438
439 ! Local Variabes
440 integer(c_int) :: nrows, ncols, p
441 logical :: ls, tb
442 type(errors) :: err
443
444 ! Initialization
445 call err%set_exit_on_error(.false.)
446 flag = la_no_error
447 tb = .false.
448 if (opb == la_transpose .or. opb == la_hermitian_transpose) tb = .true.
449 if (lside .and. tb) then
450 nrows = n
451 ncols = k
452 p = min(k, m)
453 ls = .true.
454 else if (lside .and. .not. tb) then
455 nrows = k
456 ncols = n
457 p = min(k, m)
458 ls = .true.
459 else if (.not. lside .and. tb) then
460 nrows = k
461 ncols = m
462 p = min(k, n)
463 ls = .false.
464 else
465 nrows = m
466 ncols = k
467 p = min(k, n)
468 ls = .false.
469 end if
470
471 ! Error Checking
472 if (ldb < nrows .or. ldc < m) then
473 flag = la_invalid_input_error
474 return
475 end if
476
477 ! Process
478 call diag_mtx_mult(ls, opb, alpha, a(1:p), b(1:nrows,1:ncols), &
479 beta, c(1:m,1:n))
480 if (err%has_error_occurred()) flag = err%get_error_flag()
481 end function
482
483! ------------------------------------------------------------------------------
516 function la_diag_mtx_mult_cmplx(lside, opb, m, n, k, alpha, a, b, &
517 ldb, beta, c, ldc) bind(C, name="la_diag_mtx_mult_cmplx") &
518 result(flag)
519 ! Arguments
520 logical(c_bool), intent(in), value :: lside
521 integer(c_int), intent(in), value :: opb, m, n, k, ldb, ldc
522 complex(c_double), intent(in), value :: alpha, beta
523 complex(c_double), intent(in) :: a(*), b(ldb,*)
524 complex(c_double), intent(inout) :: c(ldc,*)
525 integer(c_int) :: flag
526
527 ! Local Variabes
528 integer(c_int) :: nrows, ncols, p
529 logical :: ls, tb
530 type(errors) :: err
531
532 ! Initialization
533 call err%set_exit_on_error(.false.)
534 flag = la_no_error
535 tb = .false.
536 if (opb == la_transpose .or. opb == la_hermitian_transpose) tb = .true.
537 if (lside .and. tb) then
538 nrows = n
539 ncols = k
540 p = min(k, m)
541 ls = .true.
542 else if (lside .and. .not. tb) then
543 nrows = k
544 ncols = n
545 p = min(k, m)
546 ls = .true.
547 else if (.not. lside .and. tb) then
548 nrows = k
549 ncols = m
550 p = min(k, n)
551 ls = .false.
552 else
553 nrows = m
554 ncols = k
555 p = min(k, n)
556 ls = .false.
557 end if
558
559 ! Error Checking
560 if (ldb < nrows .or. ldc < m) then
561 flag = la_invalid_input_error
562 return
563 end if
564
565 ! Process
566 call diag_mtx_mult(ls, opb, alpha, a(1:p), b(1:nrows,1:ncols), &
567 beta, c(1:m,1:n))
568 if (err%has_error_occurred()) flag = err%get_error_flag()
569 end function
570
571! ------------------------------------------------------------------------------
588 function la_rank(m, n, a, lda, rnk) bind(C, name="la_rank") result(flag)
589 ! Arguments
590 integer(c_int), intent(in), value :: m, n, lda
591 real(c_double), intent(inout) :: a(lda,*)
592 integer(c_int), intent(out) :: rnk
593 integer(c_int) :: flag
594
595 ! Local Variables
596 type(errors) :: err
597
598 ! Input Check
599 call err%set_exit_on_error(.false.)
600 flag = la_no_error
601 if (lda < m) then
602 flag = la_invalid_input_error
603 return
604 end if
605
606 ! Process
607 rnk = mtx_rank(a(1:m,1:n), err =err)
608 if (err%has_error_occurred()) flag = err%get_error_flag()
609 end function
610
611! ------------------------------------------------------------------------------
628 function la_rank_cmplx(m, n, a, lda, rnk) bind(C, name="la_rank_cmplx") &
629 result(flag)
630 ! Arguments
631 integer(c_int), intent(in), value :: m, n, lda
632 complex(c_double), intent(inout) :: a(lda,*)
633 integer(c_int), intent(out) :: rnk
634 integer(c_int) :: flag
635
636 ! Local Variables
637 type(errors) :: err
638
639 ! Input Check
640 call err%set_exit_on_error(.false.)
641 flag = la_no_error
642 if (lda < m) then
643 flag = la_invalid_input_error
644 return
645 end if
646
647 ! Process
648 rnk = mtx_rank(a(1:m,1:n), err = err)
649 if (err%has_error_occurred()) flag = err%get_error_flag()
650 end function
651
652! ------------------------------------------------------------------------------
666 function la_det(n, a, lda, d) bind(C, name="la_det") result(flag)
667 ! Arguments
668 integer(c_int), intent(in), value :: n, lda
669 real(c_double), intent(inout) :: a(lda,*)
670 real(c_double), intent(out) :: d
671 integer(c_int) :: flag
672
673 ! Local Variables
674 type(errors) :: err
675
676 ! Error Checking
677 call err%set_exit_on_error(.false.)
678 flag = la_no_error
679 if (lda < n) then
680 flag = la_invalid_input_error
681 return
682 end if
683
684 ! Process
685 d = det(a(1:n,1:n), err = err)
686 if (err%has_error_occurred()) flag = err%get_error_flag()
687 end function
688
689! ------------------------------------------------------------------------------
703 function la_det_cmplx(n, a, lda, d) bind(C, name="la_det_cmplx") result(flag)
704 ! Arguments
705 integer(c_int), intent(in), value :: n, lda
706 complex(c_double), intent(inout) :: a(lda,*)
707 complex(c_double), intent(out) :: d
708 integer(c_int) :: flag
709
710 ! Local Variables
711 type(errors) :: err
712
713 ! Error Checking
714 call err%set_exit_on_error(.false.)
715 flag = la_no_error
716 if (lda < n) then
717 flag = la_invalid_input_error
718 return
719 end if
720
721 ! Process
722 d = det(a(1:n,1:n), err = err)
723 if (err%has_error_occurred()) flag = err%get_error_flag()
724 end function
725
726! ------------------------------------------------------------------------------
750 function la_tri_mtx_mult(upper, alpha, n, a, lda, beta, b, ldb) &
751 bind(C, name = "la_tri_mtx_mult") result(flag)
752 ! Arguments
753 logical(c_bool), intent(in), value :: upper
754 integer(c_int), intent(in), value :: n, lda, ldb
755 real(c_double), intent(in), value :: alpha, beta
756 real(c_double), intent(in) :: a(lda,*)
757 real(c_double), intent(inout) :: b(ldb,*)
758 integer(c_int) :: flag
759
760 ! Error Checking
761 flag = la_no_error
762 if (lda < n .or. ldb < n) then
763 flag = la_invalid_input_error
764 return
765 end if
766
767 ! Process
768 call tri_mtx_mult(logical(upper), alpha, a(1:n,1:n), beta, b(1:n,1:n))
769 end function
770
771! ------------------------------------------------------------------------------
795 function la_tri_mtx_mult_cmplx(upper, alpha, n, a, lda, beta, b, ldb) &
796 bind(C, name = "la_tri_mtx_mult_cmplx") result(flag)
797 ! Arguments
798 logical(c_bool), intent(in), value :: upper
799 integer(c_int), intent(in), value :: n, lda, ldb
800 complex(c_double), intent(in), value :: alpha, beta
801 complex(c_double), intent(in) :: a(lda,*)
802 complex(c_double), intent(inout) :: b(ldb,*)
803 integer(c_int) :: flag
804
805 ! Error Checking
806 flag = la_no_error
807 if (lda < n .or. ldb < n) then
808 flag = la_invalid_input_error
809 return
810 end if
811
812 ! Process
813 call tri_mtx_mult(logical(upper), alpha, a(1:n,1:n), beta, b(1:n,1:n))
814 end function
815
816! ------------------------------------------------------------------------------
834 function la_lu_factor(m, n, a, lda, ipvt) bind(C, name = "la_lu_factor") &
835 result(flag)
836 ! Arguments
837 integer(c_int), intent(in), value :: m, n, lda
838 real(c_double), intent(inout) :: a(lda,*)
839 integer(c_int), intent(out) :: ipvt(*)
840 integer(c_int) :: flag
841
842 ! Local Variables
843 type(errors) :: err
844 integer(c_int) :: mn
845
846 ! Error Checking
847 call err%set_exit_on_error(.false.)
848 flag = la_no_error
849 if (lda < m) then
850 flag = la_invalid_input_error
851 return
852 end if
853
854 ! Process
855 mn = min(m, n)
856 call lu_factor(a(1:m,1:n), ipvt(1:mn), err)
857 if (err%has_error_occurred()) then
858 flag = err%get_error_flag()
859 return
860 end if
861 end function
862
863! ------------------------------------------------------------------------------
881 function la_lu_factor_cmplx(m, n, a, lda, ipvt) &
882 bind(C, name = "la_lu_factor_cmplx") result(flag)
883 ! Arguments
884 integer(c_int), intent(in), value :: m, n, lda
885 complex(c_double), intent(inout) :: a(lda,*)
886 integer(c_int), intent(out) :: ipvt(*)
887 integer(c_int) :: flag
888
889 ! Local Variables
890 type(errors) :: err
891 integer(c_int) :: mn
892
893 ! Error Checking
894 call err%set_exit_on_error(.false.)
895 flag = la_no_error
896 if (lda < m) then
897 flag = la_invalid_input_error
898 return
899 end if
900
901 ! Process
902 mn = min(m, n)
903 call lu_factor(a(1:m,1:n), ipvt(1:mn), err)
904 if (err%has_error_occurred()) then
905 flag = err%get_error_flag()
906 return
907 end if
908 end function
909
910! ------------------------------------------------------------------------------
930 function la_form_lu(n, a, lda, ipvt, u, ldu, p, ldp) &
931 bind(C, name = "la_form_lu") result(flag)
932 ! Arguments
933 integer(c_int), intent(in), value :: n, lda, ldu, ldp
934 real(c_double), intent(inout) :: a(lda,*)
935 real(c_double), intent(out) :: u(ldu,*), p(ldp,*)
936 integer(c_int), intent(in) :: ipvt(*)
937 integer(c_int) :: flag
938
939 ! Input Checking
940 flag = la_no_error
941 if (lda < n .or. ldu < n .or. ldp < n) then
942 flag = la_invalid_input_error
943 return
944 end if
945
946 ! Process
947 call form_lu(a(1:n,1:n), ipvt(1:n), u(1:n,1:n), p(1:n,1:n))
948 end function
949
950! ------------------------------------------------------------------------------
970 function la_form_lu_cmplx(n, a, lda, ipvt, u, ldu, p, ldp) &
971 bind(C, name = "la_form_lu_cmplx") result(flag)
972 ! Arguments
973 integer(c_int), intent(in), value :: n, lda, ldu, ldp
974 complex(c_double), intent(inout) :: a(lda,*)
975 complex(c_double), intent(out) :: u(ldu,*)
976 real(c_double), intent(out) :: p(ldp,*)
977 integer(c_int), intent(in) :: ipvt(*)
978 integer(c_int) :: flag
979
980 ! Input Checking
981 flag = la_no_error
982 if (lda < n .or. ldu < n .or. ldp < n) then
983 flag = la_invalid_input_error
984 return
985 end if
986
987 ! Process
988 call form_lu(a(1:n,1:n), ipvt(1:n), u(1:n,1:n), p(1:n,1:n))
989 end function
990
991! ------------------------------------------------------------------------------
1011 function la_qr_factor(m, n, a, lda, tau) bind(C, name = "la_qr_factor") &
1012 result(flag)
1013 ! Arguments
1014 integer(c_int), intent(in), value :: m, n, lda
1015 real(c_double), intent(inout) :: a(lda,*)
1016 real(c_double), intent(out) :: tau(*)
1017 integer(c_int) :: flag
1018
1019 ! Local Variables
1020 type(errors) :: err
1021 integer(c_int) :: mn
1022
1023 ! Error Checking
1024 call err%set_exit_on_error(.false.)
1025 flag = la_no_error
1026 if (lda < m) then
1027 flag = la_invalid_input_error
1028 return
1029 end if
1030
1031 ! Process
1032 mn = min(m, n)
1033 call qr_factor(a(1:m,1:n), tau(1:mn), err = err)
1034 if (err%has_error_occurred()) then
1035 flag = err%get_error_flag()
1036 return
1037 end if
1038 end function
1039
1040! ------------------------------------------------------------------------------
1060 function la_qr_factor_cmplx(m, n, a, lda, tau) &
1061 bind(C, name = "la_qr_factor_cmplx") result(flag)
1062 ! Arguments
1063 integer(c_int), intent(in), value :: m, n, lda
1064 complex(c_double), intent(inout) :: a(lda,*)
1065 complex(c_double), intent(out) :: tau(*)
1066 integer(c_int) :: flag
1067
1068 ! Local Variables
1069 type(errors) :: err
1070 integer(c_int) :: mn
1071
1072 ! Error Checking
1073 call err%set_exit_on_error(.false.)
1074 flag = la_no_error
1075 if (lda < m) then
1076 flag = la_invalid_input_error
1077 return
1078 end if
1079
1080 ! Process
1081 mn = min(m, n)
1082 call qr_factor(a(1:m,1:n), tau(1:mn), err = err)
1083 if (err%has_error_occurred()) then
1084 flag = err%get_error_flag()
1085 return
1086 end if
1087 end function
1088
1089! ------------------------------------------------------------------------------
1113 function la_qr_factor_pvt(m, n, a, lda, tau, jpvt) &
1114 bind(C, name = "la_qr_factor_pvt") result(flag)
1115 ! Arguments
1116 integer(c_int), intent(in), value :: m, n, lda
1117 real(c_double), intent(inout) :: a(lda,*)
1118 real(c_double), intent(out) :: tau(*)
1119 integer(c_int), intent(inout) :: jpvt(*)
1120 integer(c_int) :: flag
1121
1122 ! Local Variables
1123 type(errors) :: err
1124 integer(c_int) :: mn
1125
1126 ! Error Checking
1127 call err%set_exit_on_error(.false.)
1128 flag = la_no_error
1129 if (lda < m) then
1130 flag = la_invalid_input_error
1131 return
1132 end if
1133
1134 ! Process
1135 mn = min(m, n)
1136 call qr_factor(a(1:m,1:n), tau(1:mn), jpvt(1:n), err = err)
1137 if (err%has_error_occurred()) then
1138 flag = err%get_error_flag()
1139 return
1140 end if
1141 end function
1142
1143! ------------------------------------------------------------------------------
1167 function la_qr_factor_cmplx_pvt(m, n, a, lda, tau, jpvt) &
1168 bind(C, name = "la_qr_factor_cmplx_pvt") result(flag)
1169 ! Arguments
1170 integer(c_int), intent(in), value :: m, n, lda
1171 complex(c_double), intent(inout) :: a(lda,*)
1172 complex(c_double), intent(out) :: tau(*)
1173 integer(c_int), intent(inout) :: jpvt(*)
1174 integer(c_int) :: flag
1175
1176 ! Local Variables
1177 type(errors) :: err
1178 integer(c_int) :: mn
1179
1180 ! Error Checking
1181 call err%set_exit_on_error(.false.)
1182 flag = la_no_error
1183 if (lda < m) then
1184 flag = la_invalid_input_error
1185 return
1186 end if
1187
1188 ! Process
1189 mn = min(m, n)
1190 call qr_factor(a(1:m,1:n), tau(1:mn), jpvt(1:n), err = err)
1191 if (err%has_error_occurred()) then
1192 flag = err%get_error_flag()
1193 return
1194 end if
1195 end function
1196
1197! ------------------------------------------------------------------------------
1222 function la_form_qr(fullq, m, n, r, ldr, tau, q, ldq) &
1223 bind(C, name = "la_form_qr") result(flag)
1224 ! Arguments
1225 logical(c_bool), intent(in), value :: fullq
1226 integer(c_int), intent(in), value :: m, n, ldr, ldq
1227 real(c_double), intent(inout) :: r(ldr,*)
1228 real(c_double), intent(in) :: tau(*)
1229 real(c_double), intent(out) :: q(ldq,*)
1230 integer(c_int) :: flag
1231
1232 ! Local Variables
1233 type(errors) :: err
1234 integer(c_int) :: mn, nq
1235
1236 ! Error Checking
1237 call err%set_exit_on_error(.false.)
1238 flag = la_no_error
1239 if (ldr < m .or. ldq < m) then
1240 flag = la_invalid_input_error
1241 return
1242 end if
1243
1244 ! Process
1245 mn = min(m, n)
1246 nq = m
1247 if (.not.fullq) nq = n
1248 call form_qr(r(1:m,1:n), tau(1:mn), q(1:m,1:nq), err = err)
1249 if (err%has_error_occurred()) then
1250 flag = err%get_error_flag()
1251 return
1252 end if
1253 end function
1254
1255! ------------------------------------------------------------------------------
1280 function la_form_qr_cmplx(fullq, m, n, r, ldr, tau, q, ldq) &
1281 bind(C, name = "la_form_qr_cmplx") result(flag)
1282 ! Arguments
1283 logical(c_bool), intent(in), value :: fullq
1284 integer(c_int), intent(in), value :: m, n, ldr, ldq
1285 complex(c_double), intent(inout) :: r(ldr,*)
1286 complex(c_double), intent(in) :: tau(*)
1287 complex(c_double), intent(out) :: q(ldq,*)
1288 integer(c_int) :: flag
1289
1290 ! Local Variables
1291 type(errors) :: err
1292 integer(c_int) :: mn, nq
1293
1294 ! Error Checking
1295 call err%set_exit_on_error(.false.)
1296 flag = la_no_error
1297 if (ldr < m .or. ldq < m) then
1298 flag = la_invalid_input_error
1299 return
1300 end if
1301
1302 ! Process
1303 mn = min(m, n)
1304 nq = m
1305 if (.not.fullq) nq = n
1306 call form_qr(r(1:m,1:n), tau(1:mn), q(1:m,1:nq), err = err)
1307 if (err%has_error_occurred()) then
1308 flag = err%get_error_flag()
1309 return
1310 end if
1311 end function
1312
1313! ------------------------------------------------------------------------------
1344 function la_form_qr_pvt(fullq, m, n, r, ldr, tau, pvt, q, ldq, p, ldp) &
1345 bind(C, name = "la_form_qr_pvt") result(flag)
1346 ! Arguments
1347 logical(c_bool), intent(in), value :: fullq
1348 integer(c_int), intent(in), value :: m, n, ldr, ldq, ldp
1349 real(c_double), intent(inout) :: r(ldr,*)
1350 real(c_double), intent(in) :: tau(*)
1351 integer(c_int), intent(in) :: pvt(*)
1352 real(c_double), intent(out) :: q(ldq,*), p(ldp,*)
1353 integer(c_int) :: flag
1354
1355 ! Local Variables
1356 type(errors) :: err
1357 integer(c_int) :: mn, nq
1358
1359 ! Error Checking
1360 call err%set_exit_on_error(.false.)
1361 flag = la_no_error
1362 if (ldr < m .or. ldq < m .or. ldp < n) then
1363 flag = la_invalid_input_error
1364 return
1365 end if
1366
1367 ! Process
1368 mn = min(m, n)
1369 nq = m
1370 if (.not.fullq) nq = n
1371 call form_qr(r(1:m,1:n), tau(1:mn), pvt(1:n), q(1:m,1:nq), p(1:n,1:n), &
1372 err = err)
1373 if (err%has_error_occurred()) then
1374 flag = err%get_error_flag()
1375 return
1376 end if
1377 end function
1378
1379! ------------------------------------------------------------------------------
1410 function la_form_qr_cmplx_pvt(fullq, m, n, r, ldr, tau, pvt, q, ldq, p, &
1411 ldp) bind(C, name = "la_form_qr_cmplx_pvt") result(flag)
1412 ! Arguments
1413 logical(c_bool), intent(in), value :: fullq
1414 integer(c_int), intent(in), value :: m, n, ldr, ldq, ldp
1415 complex(c_double), intent(inout) :: r(ldr,*)
1416 complex(c_double), intent(in) :: tau(*)
1417 integer(c_int), intent(in) :: pvt(*)
1418 complex(c_double), intent(out) :: q(ldq,*), p(ldp,*)
1419 integer(c_int) :: flag
1420
1421 ! Local Variables
1422 type(errors) :: err
1423 integer(c_int) :: mn, nq
1424
1425 ! Error Checking
1426 call err%set_exit_on_error(.false.)
1427 flag = la_no_error
1428 if (ldr < m .or. ldq < m .or. ldp < n) then
1429 flag = la_invalid_input_error
1430 return
1431 end if
1432
1433 ! Process
1434 mn = min(m, n)
1435 nq = m
1436 if (.not.fullq) nq = n
1437 call form_qr(r(1:m,1:n), tau(1:mn), pvt(1:n), q(1:m,1:nq), p(1:n,1:n), &
1438 err = err)
1439 if (err%has_error_occurred()) then
1440 flag = err%get_error_flag()
1441 return
1442 end if
1443 end function
1444
1445! ------------------------------------------------------------------------------
1473 function la_mult_qr(lside, trans, m, n, k, a, lda, tau, c, ldc) &
1474 bind(C, name = "la_mult_qr") result(flag)
1475 ! Local Variables
1476 logical(c_bool), intent(in), value :: lside, trans
1477 integer(c_int), intent(in), value :: m, n, k, lda, ldc
1478 real(c_double), intent(inout) :: a(lda,*), c(ldc,*)
1479 real(c_double), intent(in) :: tau(*)
1480 integer(c_int) :: flag
1481
1482 ! Local Variables
1483 type(errors) :: err
1484 integer(c_int) :: ma, na
1485
1486 ! Initialization
1487 if (lside) then
1488 ma = m
1489 na = m
1490 else
1491 ma = n
1492 na = n
1493 end if
1494
1495 ! Error Checking
1496 call err%set_exit_on_error(.false.)
1497 flag = la_no_error
1498 if (lda < ma .or. ldc < m) then
1499 flag = la_invalid_input_error
1500 return
1501 end if
1502 if (k > na .or. k < 0) then
1503 flag = la_invalid_input_error
1504 return
1505 end if
1506
1507 ! Process
1508 call mult_qr(logical(lside), logical(trans), a(1:ma,1:k), tau(1:k), &
1509 c(1:m,1:n), err = err)
1510 if (err%has_error_occurred()) then
1511 flag = err%get_error_flag()
1512 return
1513 end if
1514 end function
1515
1516! ------------------------------------------------------------------------------
1544 function la_mult_qr_cmplx(lside, trans, m, n, k, a, lda, tau, c, ldc) &
1545 bind(C, name = "la_mult_qr_cmplx") result(flag)
1546 ! Local Variables
1547 logical(c_bool), intent(in), value :: lside, trans
1548 integer(c_int), intent(in), value :: m, n, k, lda, ldc
1549 complex(c_double), intent(inout) :: a(lda,*), c(ldc,*)
1550 complex(c_double), intent(in) :: tau(*)
1551 integer(c_int) :: flag
1552
1553 ! Local Variables
1554 type(errors) :: err
1555 integer(c_int) :: ma, na
1556
1557 ! Initialization
1558 if (lside) then
1559 ma = m
1560 na = m
1561 else
1562 ma = n
1563 na = n
1564 end if
1565
1566 ! Error Checking
1567 call err%set_exit_on_error(.false.)
1568 flag = la_no_error
1569 if (lda < ma .or. ldc < m) then
1570 flag = la_invalid_input_error
1571 return
1572 end if
1573 if (k > na .or. k < 0) then
1574 flag = la_invalid_input_error
1575 return
1576 end if
1577
1578 ! Process
1579 call mult_qr(logical(lside), logical(trans), a(1:ma,1:k), tau(1:k), &
1580 c(1:m,1:n), err = err)
1581 if (err%has_error_occurred()) then
1582 flag = err%get_error_flag()
1583 return
1584 end if
1585 end function
1586
1587! ------------------------------------------------------------------------------
1609 function la_qr_rank1_update(m, n, q, ldq, r, ldr, u, v) &
1610 bind(C, name = "la_qr_rank1_update") result(flag)
1611 ! Arguments
1612 integer(c_int), intent(in), value :: m, n, ldq, ldr
1613 real(c_double), intent(inout) :: q(ldq,*), r(ldr,*), u(*), v(*)
1614 integer(c_int) :: flag
1615
1616 ! Local Variables
1617 type(errors) :: err
1618 integer(c_int) :: mn
1619
1620 ! Error Checking
1621 call err%set_exit_on_error(.false.)
1622 flag = la_no_error
1623 if (ldq < m .or. ldr < m) then
1624 flag = la_invalid_input_error
1625 return
1626 end if
1627
1628 ! Process
1629 mn = min(m, n)
1630 call qr_rank1_update(q(1:m,1:mn), r(1:m,1:n), u(1:m), v(1:n), err = err)
1631 if (err%has_error_occurred()) then
1632 flag = err%get_error_flag()
1633 return
1634 end if
1635 end function
1636
1637! ------------------------------------------------------------------------------
1659 function la_qr_rank1_update_cmplx(m, n, q, ldq, r, ldr, u, v) &
1660 bind(C, name = "la_qr_rank1_update_cmplx") result(flag)
1661 ! Arguments
1662 integer(c_int), intent(in), value :: m, n, ldq, ldr
1663 complex(c_double), intent(inout) :: q(ldq,*), r(ldr,*), u(*), v(*)
1664 integer(c_int) :: flag
1665
1666 ! Local Variables
1667 type(errors) :: err
1668 integer(c_int) :: mn
1669
1670 ! Error Checking
1671 call err%set_exit_on_error(.false.)
1672 flag = la_no_error
1673 if (ldq < m .or. ldr < m) then
1674 flag = la_invalid_input_error
1675 return
1676 end if
1677
1678 ! Process
1679 mn = min(m, n)
1680 call qr_rank1_update(q(1:m,1:mn), r(1:m,1:n), u(1:m), v(1:n), err = err)
1681 if (err%has_error_occurred()) then
1682 flag = err%get_error_flag()
1683 return
1684 end if
1685 end function
1686
1687! ------------------------------------------------------------------------------
1704 function la_cholesky_factor(upper, n, a, lda) &
1705 bind(C, name = "la_cholesky_factor") result(flag)
1706 ! Arguments
1707 logical(c_bool), intent(in), value :: upper
1708 integer(c_int), intent(in), value :: n, lda
1709 real(c_double), intent(inout) :: a(lda,*)
1710 integer(c_int) :: flag
1711
1712 ! Local Variables
1713 type(errors) :: err
1714
1715 ! Error Checking
1716 call err%set_exit_on_error(.false.)
1717 flag = la_no_error
1718 if (lda < n) then
1719 flag = la_invalid_input_error
1720 return
1721 end if
1722
1723 ! Process
1724 call cholesky_factor(a(1:n,1:n), logical(upper), err = err)
1725 if (err%has_error_occurred()) then
1726 flag = err%get_error_flag()
1727 return
1728 end if
1729 end function
1730
1731! ------------------------------------------------------------------------------
1748 function la_cholesky_factor_cmplx(upper, n, a, lda) &
1749 bind(C, name = "la_cholesky_factor_cmplx") result(flag)
1750 ! Arguments
1751 logical(c_bool), intent(in), value :: upper
1752 integer(c_int), intent(in), value :: n, lda
1753 complex(c_double), intent(inout) :: a(lda,*)
1754 integer(c_int) :: flag
1755
1756 ! Local Variables
1757 type(errors) :: err
1758
1759 ! Error Checking
1760 call err%set_exit_on_error(.false.)
1761 flag = la_no_error
1762 if (lda < n) then
1763 flag = la_invalid_input_error
1764 return
1765 end if
1766
1767 ! Process
1768 call cholesky_factor(a(1:n,1:n), logical(upper), err = err)
1769 if (err%has_error_occurred()) then
1770 flag = err%get_error_flag()
1771 return
1772 end if
1773 end function
1774
1775! ------------------------------------------------------------------------------
1791 function la_cholesky_rank1_update(n, r, ldr, u) &
1792 bind(C, name = "la_cholesky_rank1_update") result(flag)
1793 ! Arguments
1794 integer(c_int), intent(in), value :: n, ldr
1795 real(c_double), intent(inout) :: r(ldr,*), u(*)
1796 integer(c_int) :: flag
1797
1798 ! Local Variables
1799 type(errors) :: err
1800
1801 ! Error Checking
1802 call err%set_exit_on_error(.false.)
1803 flag = la_no_error
1804 if (ldr < n) then
1805 flag = la_invalid_input_error
1806 return
1807 end if
1808
1809 ! Process
1810 call cholesky_rank1_update(r(1:n,1:n), u(1:n), err = err)
1811 if (err%has_error_occurred()) then
1812 flag = err%get_error_flag()
1813 return
1814 end if
1815 end function
1816
1817! ------------------------------------------------------------------------------
1833 function la_cholesky_rank1_update_cmplx(n, r, ldr, u) &
1834 bind(C, name = "la_cholesky_rank1_update_cmplx") result(flag)
1835 ! Arguments
1836 integer(c_int), intent(in), value :: n, ldr
1837 complex(c_double), intent(inout) :: r(ldr,*), u(*)
1838 integer(c_int) :: flag
1839
1840 ! Local Variables
1841 type(errors) :: err
1842
1843 ! Error Checking
1844 call err%set_exit_on_error(.false.)
1845 flag = la_no_error
1846 if (ldr < n) then
1847 flag = la_invalid_input_error
1848 return
1849 end if
1850
1851 ! Process
1852 call cholesky_rank1_update(r(1:n,1:n), u(1:n), err = err)
1853 if (err%has_error_occurred()) then
1854 flag = err%get_error_flag()
1855 return
1856 end if
1857 end function
1858
1859! ------------------------------------------------------------------------------
1877 function la_cholesky_rank1_downdate(n, r, ldr, u) &
1878 bind(C, name = "la_cholesky_rank1_downdate") result(flag)
1879 ! Arguments
1880 integer(c_int), intent(in), value :: n, ldr
1881 real(c_double), intent(inout) :: r(ldr,*), u(*)
1882 integer(c_int) :: flag
1883
1884 ! Local Variables
1885 type(errors) :: err
1886
1887 ! Error Checking
1888 call err%set_exit_on_error(.false.)
1889 flag = la_no_error
1890 if (ldr < n) then
1891 flag = la_invalid_input_error
1892 return
1893 end if
1894
1895 ! Process
1896 call cholesky_rank1_downdate(r(1:n,1:n), u(1:n), err = err)
1897 if (err%has_error_occurred()) then
1898 flag = err%get_error_flag()
1899 return
1900 end if
1901 end function
1902
1903! ------------------------------------------------------------------------------
1921 function la_cholesky_rank1_downdate_cmplx(n, r, ldr, u) &
1922 bind(C, name = "la_cholesky_rank1_downdate_cmplx") result(flag)
1923 ! Arguments
1924 integer(c_int), intent(in), value :: n, ldr
1925 complex(c_double), intent(inout) :: r(ldr,*), u(*)
1926 integer(c_int) :: flag
1927
1928 ! Local Variables
1929 type(errors) :: err
1930
1931 ! Error Checking
1932 call err%set_exit_on_error(.false.)
1933 flag = la_no_error
1934 if (ldr < n) then
1935 flag = la_invalid_input_error
1936 return
1937 end if
1938
1939 ! Process
1940 call cholesky_rank1_downdate(r(1:n,1:n), u(1:n), err = err)
1941 if (err%has_error_occurred()) then
1942 flag = err%get_error_flag()
1943 return
1944 end if
1945 end function
1946
1947! ------------------------------------------------------------------------------
1975 function la_svd(m, n, a, lda, s, u, ldu, vt, ldv) &
1976 bind(C, name = "la_svd") result(flag)
1977 ! Arguments
1978 integer(c_int), intent(in), value :: m, n, lda, ldu, ldv
1979 real(c_double), intent(inout) :: a(lda,*)
1980 real(c_double), intent(out) :: s(*), u(ldu,*), vt(ldv,*)
1981 integer(c_int) :: flag
1982
1983 ! Local Variables
1984 type(errors) :: err
1985 integer(c_int) :: mn
1986
1987 ! Error Checking
1988 call err%set_exit_on_error(.false.)
1989 flag = la_no_error
1990 if (lda < m .or. ldu < m .or. ldv < n) then
1991 flag = la_invalid_input_error
1992 return
1993 end if
1994
1995 ! Process
1996 mn = min(m, n)
1997 call svd(a(1:m,1:n), s(1:mn), u(1:m,1:m), vt(1:n,1:n), err = err)
1998 if (err%has_error_occurred()) then
1999 flag = err%get_error_flag()
2000 return
2001 end if
2002 end function
2003
2004! ------------------------------------------------------------------------------
2032 function la_svd_cmplx(m, n, a, lda, s, u, ldu, vt, ldv) &
2033 bind(C, name = "la_svd_cmplx") result(flag)
2034 ! Arguments
2035 integer(c_int), intent(in), value :: m, n, lda, ldu, ldv
2036 complex(c_double), intent(inout) :: a(lda,*)
2037 real(c_double), intent(out) :: s(*)
2038 complex(c_double), intent(out) :: u(ldu,*), vt(ldv,*)
2039 integer(c_int) :: flag
2040
2041 ! Local Variables
2042 type(errors) :: err
2043 integer(c_int) :: mn
2044
2045 ! Error Checking
2046 call err%set_exit_on_error(.false.)
2047 flag = la_no_error
2048 if (lda < m .or. ldu < m .or. ldv < n) then
2049 flag = la_invalid_input_error
2050 return
2051 end if
2052
2053 ! Process
2054 mn = min(m, n)
2055 call svd(a(1:m,1:n), s(1:mn), u(1:m,1:m), vt(1:n,1:n), err = err)
2056 if (err%has_error_occurred()) then
2057 flag = err%get_error_flag()
2058 return
2059 end if
2060 end function
2061
2062! ------------------------------------------------------------------------------
2089 function la_solve_tri_mtx(lside, upper, trans, nounit, m, n, alpha, a, &
2090 lda, b, ldb) bind(C, name = "la_solve_tri_mtx") result(flag)
2091 ! Arguments
2092 logical(c_bool), intent(in), value :: lside, upper, trans, nounit
2093 integer(c_int), intent(in), value :: m, n, lda, ldb
2094 real(c_double), intent(in), value :: alpha
2095 real(c_double), intent(in) :: a(lda,*)
2096 real(c_double), intent(inout) :: b(ldb,*)
2097 integer(c_int) :: flag
2098
2099 ! Local Variables
2100 type(errors) :: err
2101 integer(c_int) :: ma
2102
2103 ! Initialization
2104 if (lside) then
2105 ma = m
2106 else
2107 ma = n
2108 end if
2109
2110 ! Error Checking
2111 call err%set_exit_on_error(.false.)
2112 flag = la_no_error
2113 if (lda < ma .or. ldb < m) then
2114 flag = la_invalid_input_error
2115 return
2116 end if
2117
2118 ! Process
2119 call solve_triangular_system(logical(lside), logical(upper), &
2120 logical(trans), logical(nounit), alpha, a(1:ma,1:ma), b(1:m,1:n))
2121 end function
2122
2123! ------------------------------------------------------------------------------
2150 function la_solve_tri_mtx_cmplx(lside, upper, trans, nounit, m, n, &
2151 alpha, a, lda, b, ldb) &
2152 bind(C, name = "la_solve_tri_mtx_cmplx") result(flag)
2153 ! Arguments
2154 logical(c_bool), intent(in), value :: lside, upper, trans, nounit
2155 integer(c_int), intent(in), value :: m, n, lda, ldb
2156 complex(c_double), intent(in), value :: alpha
2157 complex(c_double), intent(in) :: a(lda,*)
2158 complex(c_double), intent(inout) :: b(ldb,*)
2159 integer(c_int) :: flag
2160
2161 ! Local Variables
2162 type(errors) :: err
2163 integer(c_int) :: ma
2164
2165 ! Initialization
2166 if (lside) then
2167 ma = m
2168 else
2169 ma = n
2170 end if
2171
2172 ! Error Checking
2173 call err%set_exit_on_error(.false.)
2174 flag = la_no_error
2175 if (lda < ma .or. ldb < m) then
2176 flag = la_invalid_input_error
2177 return
2178 end if
2179
2180 ! Process
2181 call solve_triangular_system(logical(lside), logical(upper), &
2182 logical(trans), logical(nounit), alpha, a(1:ma,1:ma), b(1:m,1:n))
2183 end function
2184
2185! ------------------------------------------------------------------------------
2200 function la_solve_lu(m, n, a, lda, ipvt, b, ldb) &
2201 bind(C, name = "la_solve_lu") result(flag)
2202 ! Arguments
2203 integer(c_int), intent(in), value :: m, n, lda, ldb
2204 real(c_double), intent(in) :: a(lda,*)
2205 integer(c_int), intent(in) :: ipvt(*)
2206 real(c_double), intent(inout) :: b(ldb,*)
2207 integer(c_int) :: flag
2208
2209 ! Local Variables
2210 type(errors) :: err
2211
2212 ! Error Checking
2213 call err%set_exit_on_error(.false.)
2214 flag = la_no_error
2215 if (lda < m .or. ldb < m) then
2216 flag = la_invalid_input_error
2217 return
2218 end if
2219
2220 ! Process
2221 call solve_lu(a(1:m,1:m), ipvt(1:m), b(1:m,1:n))
2222 end function
2223
2224! ------------------------------------------------------------------------------
2239 function la_solve_lu_cmplx(m, n, a, lda, ipvt, b, ldb) &
2240 bind(C, name = "la_solve_lu_cmplx") result(flag)
2241 ! Arguments
2242 integer(c_int), intent(in), value :: m, n, lda, ldb
2243 complex(c_double), intent(in) :: a(lda,*)
2244 integer(c_int), intent(in) :: ipvt(*)
2245 complex(c_double), intent(inout) :: b(ldb,*)
2246 integer(c_int) :: flag
2247
2248 ! Local Variables
2249 type(errors) :: err
2250
2251 ! Error Checking
2252 call err%set_exit_on_error(.false.)
2253 flag = la_no_error
2254 if (lda < m .or. ldb < m) then
2255 flag = la_invalid_input_error
2256 return
2257 end if
2258
2259 ! Process
2260 call solve_lu(a(1:m,1:m), ipvt(1:m), b(1:m,1:n))
2261 end function
2262
2263! ------------------------------------------------------------------------------
2285 function la_solve_qr(m, n, k, a, lda, tau, b, ldb) &
2286 bind(C, name = "la_solve_qr") result(flag)
2287 ! Arguments
2288 integer(c_int), intent(in), value :: m, n, k, lda, ldb
2289 real(c_double), intent(inout) :: a(lda,*), b(ldb,*)
2290 real(c_double), intent(in) :: tau(*)
2291 integer(c_int) :: flag
2292
2293 ! Local Variables
2294 type(errors) :: err
2295 integer(c_int) :: minmn
2296
2297 ! Error Checking
2298 call err%set_exit_on_error(.false.)
2299 flag = la_no_error
2300 if (lda < m .or. ldb < m .or. m < n) then
2301 flag = la_invalid_input_error
2302 return
2303 end if
2304
2305 ! Process
2306 minmn = min(m, n)
2307 call solve_qr(a(1:m,1:n), tau(1:minmn), b(1:m,1:k), err = err)
2308 if (err%has_error_occurred()) then
2309 flag = err%get_error_flag()
2310 return
2311 end if
2312 end function
2313
2314! ------------------------------------------------------------------------------
2336 function la_solve_qr_cmplx(m, n, k, a, lda, tau, b, ldb) &
2337 bind(C, name = "la_solve_qr_cmplx") result(flag)
2338 ! Arguments
2339 integer(c_int), intent(in), value :: m, n, k, lda, ldb
2340 complex(c_double), intent(inout) :: a(lda,*), b(ldb,*)
2341 complex(c_double), intent(in) :: tau(*)
2342 integer(c_int) :: flag
2343
2344 ! Local Variables
2345 type(errors) :: err
2346 integer(c_int) :: minmn
2347
2348 ! Error Checking
2349 call err%set_exit_on_error(.false.)
2350 flag = la_no_error
2351 if (lda < m .or. ldb < m .or. m < n) then
2352 flag = la_invalid_input_error
2353 return
2354 end if
2355
2356 ! Process
2357 minmn = min(m, n)
2358 call solve_qr(a(1:m,1:n), tau(1:minmn), b(1:m,1:k), err = err)
2359 if (err%has_error_occurred()) then
2360 flag = err%get_error_flag()
2361 return
2362 end if
2363 end function
2364
2365! ------------------------------------------------------------------------------
2387 function la_solve_qr_pvt(m, n, k, a, lda, tau, jpvt, b, ldb) &
2388 bind(C, name = "la_solve_qr_pvt") result(flag)
2389 ! Arguments
2390 integer(c_int), intent(in), value :: m, n, k, lda, ldb
2391 real(c_double), intent(inout) :: a(lda,*), b(ldb,*)
2392 real(c_double), intent(in) :: tau(*)
2393 integer(c_int), intent(in) :: jpvt(*)
2394 integer(c_int) :: flag
2395
2396 ! Local Variables
2397 type(errors) :: err
2398 integer(c_int) :: minmn, maxmn
2399
2400 ! Error Checking
2401 minmn = min(m, n)
2402 maxmn = max(m, n)
2403 call err%set_exit_on_error(.false.)
2404 flag = la_no_error
2405 if (lda < m .or. ldb < maxmn) then
2406 flag = la_invalid_input_error
2407 return
2408 end if
2409
2410 ! Process
2411 call solve_qr(a(1:m,1:n), tau(1:minmn), jpvt(1:n), b(1:maxmn,1:k), &
2412 err = err)
2413 if (err%has_error_occurred()) then
2414 flag = err%get_error_flag()
2415 return
2416 end if
2417 end function
2418
2419! ------------------------------------------------------------------------------
2441 function la_solve_qr_cmplx_pvt(m, n, k, a, lda, tau, jpvt, b, ldb) &
2442 bind(C, name = "la_solve_qr_cmplx_pvt") result(flag)
2443 ! Arguments
2444 integer(c_int), intent(in), value :: m, n, k, lda, ldb
2445 complex(c_double), intent(inout) :: a(lda,*), b(ldb,*)
2446 complex(c_double), intent(in) :: tau(*)
2447 integer(c_int), intent(in) :: jpvt(*)
2448 integer(c_int) :: flag
2449
2450 ! Local Variables
2451 type(errors) :: err
2452 integer(c_int) :: minmn, maxmn
2453
2454 ! Error Checking
2455 minmn = min(m, n)
2456 maxmn = max(m, n)
2457 call err%set_exit_on_error(.false.)
2458 flag = la_no_error
2459 if (lda < m .or. ldb < maxmn) then
2460 flag = la_invalid_input_error
2461 return
2462 end if
2463
2464 ! Process
2465 call solve_qr(a(1:m,1:n), tau(1:minmn), jpvt(1:n), b(1:maxmn,1:k), &
2466 err = err)
2467 if (err%has_error_occurred()) then
2468 flag = err%get_error_flag()
2469 return
2470 end if
2471 end function
2472
2473! ------------------------------------------------------------------------------
2490 function la_solve_cholesky(upper, m, n, a, lda, b, ldb) &
2491 bind(C, name = "la_solve_cholesky") result(flag)
2492 ! Arguments
2493 logical(c_bool), intent(in), value :: upper
2494 integer(c_int), intent(in), value :: m, n, lda, ldb
2495 real(c_double), intent(in) :: a(lda,*)
2496 real(c_double), intent(inout) :: b(ldb,*)
2497 integer(c_int) :: flag
2498
2499 ! Local Variables
2500 type(errors) :: err
2501
2502 ! Error Checking
2503 call err%set_exit_on_error(.false.)
2504 flag = la_no_error
2505 if (lda < m .or. ldb < m) then
2506 flag = la_invalid_input_error
2507 return
2508 end if
2509
2510 ! Process
2511 call solve_cholesky(logical(upper), a(1:m,1:m), b(1:m,1:n))
2512 end function
2513
2514! ------------------------------------------------------------------------------
2531 function la_solve_cholesky_cmplx(upper, m, n, a, lda, b, ldb) &
2532 bind(C, name = "la_solve_cholesky_cmplx") result(flag)
2533 ! Arguments
2534 logical(c_bool), intent(in), value :: upper
2535 integer(c_int), intent(in), value :: m, n, lda, ldb
2536 complex(c_double), intent(in) :: a(lda,*)
2537 complex(c_double), intent(inout) :: b(ldb,*)
2538 integer(c_int) :: flag
2539
2540 ! Local Variables
2541 type(errors) :: err
2542
2543 ! Error Checking
2544 call err%set_exit_on_error(.false.)
2545 flag = la_no_error
2546 if (lda < m .or. ldb < m) then
2547 flag = la_invalid_input_error
2548 return
2549 end if
2550
2551 ! Process
2552 call solve_cholesky(logical(upper), a(1:m,1:m), b(1:m,1:n))
2553 end function
2554
2555! ------------------------------------------------------------------------------
2579 function la_solve_least_squares(m, n, k, a, lda, b, ldb) &
2580 bind(C, name = "la_solve_least_squares") result(flag)
2581 ! Arguments
2582 integer(c_int), intent(in), value :: m, n, k, lda, ldb
2583 real(c_double), intent(inout) :: a(lda,*), b(ldb,*)
2584 integer(c_int) :: flag
2585
2586 ! Local Variables
2587 type(errors) :: err
2588 integer(c_int) :: maxmn
2589
2590 ! Error Checking
2591 maxmn = max(m, n)
2592 call err%set_exit_on_error(.false.)
2593 flag = la_no_error
2594 if (lda < m .or. ldb < maxmn) then
2595 flag = la_invalid_input_error
2596 return
2597 end if
2598
2599 ! Process
2600 call solve_least_squares(a(1:m,1:n), b(1:maxmn,1:k), err = err)
2601 if (err%has_error_occurred()) then
2602 flag = err%get_error_flag()
2603 return
2604 end if
2605 end function
2606
2607! ------------------------------------------------------------------------------
2631 function la_solve_least_squares_cmplx(m, n, k, a, lda, b, ldb) &
2632 bind(C, name = "la_solve_least_squares_cmplx") result(flag)
2633 ! Arguments
2634 integer(c_int), intent(in), value :: m, n, k, lda, ldb
2635 complex(c_double), intent(inout) :: a(lda,*), b(ldb,*)
2636 integer(c_int) :: flag
2637
2638 ! Local Variables
2639 type(errors) :: err
2640 integer(c_int) :: maxmn
2641
2642 ! Error Checking
2643 maxmn = max(m, n)
2644 call err%set_exit_on_error(.false.)
2645 flag = la_no_error
2646 if (lda < m .or. ldb < maxmn) then
2647 flag = la_invalid_input_error
2648 return
2649 end if
2650
2651 ! Process
2652 call solve_least_squares(a(1:m,1:n), b(1:maxmn,1:k), err = err)
2653 if (err%has_error_occurred()) then
2654 flag = err%get_error_flag()
2655 return
2656 end if
2657 end function
2658
2659! ------------------------------------------------------------------------------
2671 function la_inverse(n, a, lda) bind(C, name = "la_inverse") result(flag)
2672 ! Arguments
2673 integer(c_int), intent(in), value :: n, lda
2674 real(c_double), intent(inout) :: a(lda,*)
2675 integer(c_int) :: flag
2676
2677 ! Local Variables
2678 type(errors) :: err
2679
2680 ! Error Checking
2681 call err%set_exit_on_error(.false.)
2682 flag = la_no_error
2683 if (lda < n) then
2684 flag = la_invalid_input_error
2685 return
2686 end if
2687
2688 ! Process
2689 call mtx_inverse(a(1:n,1:n), err = err)
2690 if (err%has_error_occurred()) then
2691 flag = err%get_error_flag()
2692 return
2693 end if
2694 end function
2695
2696! ------------------------------------------------------------------------------
2708 function la_inverse_cmplx(n, a, lda) bind(C, name = "la_inverse_cmplx") &
2709 result(flag)
2710 ! Arguments
2711 integer(c_int), intent(in), value :: n, lda
2712 complex(c_double), intent(inout) :: a(lda,*)
2713 integer(c_int) :: flag
2714
2715 ! Local Variables
2716 type(errors) :: err
2717
2718 ! Error Checking
2719 call err%set_exit_on_error(.false.)
2720 flag = la_no_error
2721 if (lda < n) then
2722 flag = la_invalid_input_error
2723 return
2724 end if
2725
2726 ! Process
2727 call mtx_inverse(a(1:n,1:n), err = err)
2728 if (err%has_error_occurred()) then
2729 flag = err%get_error_flag()
2730 return
2731 end if
2732 end function
2733
2734! ------------------------------------------------------------------------------
2750 function la_pinverse(m, n, a, lda, ainv, ldai) &
2751 bind(C, name = "la_pinverse") result(flag)
2752 ! Arguments
2753 integer(c_int), intent(in), value :: m, n, lda, ldai
2754 real(c_double), intent(inout) :: a(lda,*)
2755 real(c_double), intent(out) :: ainv(ldai,*)
2756 integer(c_int) :: flag
2757
2758 ! Local Variables
2759 type(errors) :: err
2760
2761 ! Error Checking
2762 call err%set_exit_on_error(.false.)
2763 flag = la_no_error
2764 if (lda < m .or. ldai < n) then
2765 flag = la_invalid_input_error
2766 return
2767 end if
2768
2769 ! Process
2770 call mtx_pinverse(a(1:m,1:n), ainv(1:n,1:m), err = err)
2771 if (err%has_error_occurred()) then
2772 flag = err%get_error_flag()
2773 return
2774 end if
2775 end function
2776
2777! ------------------------------------------------------------------------------
2793 function la_pinverse_cmplx(m, n, a, lda, ainv, ldai) &
2794 bind(C, name = "la_pinverse_cmplx") result(flag)
2795 ! Arguments
2796 integer(c_int), intent(in), value :: m, n, lda, ldai
2797 complex(c_double), intent(inout) :: a(lda,*)
2798 complex(c_double), intent(out) :: ainv(ldai,*)
2799 integer(c_int) :: flag
2800
2801 ! Local Variables
2802 type(errors) :: err
2803
2804 ! Error Checking
2805 call err%set_exit_on_error(.false.)
2806 flag = la_no_error
2807 if (lda < m .or. ldai < n) then
2808 flag = la_invalid_input_error
2809 return
2810 end if
2811
2812 ! Process
2813 call mtx_pinverse(a(1:m,1:n), ainv(1:n,1:m), err = err)
2814 if (err%has_error_occurred()) then
2815 flag = err%get_error_flag()
2816 return
2817 end if
2818 end function
2819
2820! ------------------------------------------------------------------------------
2842 function la_eigen_symm(vecs, n, a, lda, vals) &
2843 bind(C, name = "la_eigen_symm") result(flag)
2844 ! Arguments
2845 logical(c_bool), intent(in), value :: vecs
2846 integer(c_int), intent(in), value :: n, lda
2847 real(c_double), intent(inout) :: a(lda,*)
2848 real(c_double), intent(out) :: vals(*)
2849 integer(c_int) :: flag
2850
2851 ! Local Variables
2852 type(errors) :: err
2853
2854 ! Error Checking
2855 call err%set_exit_on_error(.false.)
2856 flag = la_no_error
2857 if (lda < n) then
2858 flag = la_invalid_input_error
2859 return
2860 end if
2861
2862 ! Process
2863 call eigen(logical(vecs), a(1:n,1:n), vals(1:n), err = err)
2864 if (err%has_error_occurred()) then
2865 flag = err%get_error_flag()
2866 return
2867 end if
2868 end function
2869
2870! ------------------------------------------------------------------------------
2891 function la_eigen_asymm(vecs, n, a, lda, vals, v, ldv) &
2892 bind(C, name = "la_eigen_asymm") result(flag)
2893 ! Arguments
2894 logical(c_bool), intent(in), value :: vecs
2895 integer(c_int), intent(in), value :: n, lda, ldv
2896 real(c_double), intent(inout) :: a(lda,*)
2897 complex(c_double), intent(out) :: vals(*), v(ldv,*)
2898 integer(c_int) :: flag
2899
2900 ! Local Variables
2901 type(errors) :: err
2902
2903 ! Error Checking
2904 call err%set_exit_on_error(.false.)
2905 flag = la_no_error
2906 if (vecs) then
2907 if (lda < n .or. ldv < n) then
2908 flag = la_invalid_input_error
2909 return
2910 end if
2911 else
2912 if (lda < n) then
2913 flag = la_invalid_input_error
2914 return
2915 end if
2916 end if
2917
2918 ! Process
2919 if (vecs) then
2920 call eigen(a(1:n,1:n), vals(1:n), v(1:n,1:n), err = err)
2921 else
2922 call eigen(a(1:n,1:n), vals(1:n))
2923 end if
2924 if (err%has_error_occurred()) then
2925 flag = err%get_error_flag()
2926 return
2927 end if
2928 end function
2929
2930! ------------------------------------------------------------------------------
2964 function la_eigen_gen(vecs, n, a, lda, b, ldb, alpha, beta, v, ldv) &
2965 bind(C, name = "la_eigen_gen") result(flag)
2966 ! Arguments
2967 logical(c_bool), intent(in), value :: vecs
2968 integer(c_int), intent(in), value :: n, lda, ldb, ldv
2969 real(c_double), intent(inout) :: a(lda,*), b(ldb,*)
2970 real(c_double), intent(out) :: beta(*)
2971 complex(c_double), intent(out) :: alpha(*), v(ldv,*)
2972 integer(c_int) :: flag
2973
2974 ! Local Variables
2975 type(errors) :: err
2976
2977 ! Error Checking
2978 call err%set_exit_on_error(.false.)
2979 flag = la_no_error
2980 if (vecs) then
2981 if (lda < n .or. ldb < n .or. ldv < n) then
2982 flag = la_invalid_input_error
2983 return
2984 end if
2985 else
2986 if (lda < n .or. ldb < n) then
2987 flag = la_invalid_input_error
2988 return
2989 end if
2990 end if
2991
2992 ! Process
2993 if (vecs) then
2994 call eigen(a(1:n,1:n), b(1:n,1:n), alpha(1:n), beta(1:n), &
2995 v(1:n,1:n), err = err)
2996 else
2997 call eigen(a(1:n,1:n), b(1:n,1:n), alpha(1:n), beta(1:n), err = err)
2998 end if
2999 if (err%has_error_occurred()) then
3000 flag = err%get_error_flag()
3001 return
3002 end if
3003 end function
3004
3005! ------------------------------------------------------------------------------
3026 function la_eigen_cmplx(vecs, n, a, lda, vals, v, ldv) &
3027 bind(C, name = "la_eigen_cmplx") result(flag)
3028 ! Arguments
3029 logical(c_bool), intent(in), value :: vecs
3030 integer(c_int), intent(in), value :: n, lda, ldv
3031 complex(c_double), intent(inout) :: a(lda,*)
3032 complex(c_double), intent(out) :: vals(*), v(ldv,*)
3033 integer(c_int) :: flag
3034
3035 ! Local Variables
3036 type(errors) :: err
3037
3038 ! Error Checking
3039 call err%set_exit_on_error(.false.)
3040 flag = la_no_error
3041 if (vecs) then
3042 if (lda < n .or. ldv < n) then
3043 flag = la_invalid_input_error
3044 return
3045 end if
3046 else
3047 if (lda < n) then
3048 flag = la_invalid_input_error
3049 return
3050 end if
3051 end if
3052
3053 ! Process
3054 if (vecs) then
3055 call eigen(a(1:n,1:n), vals(1:n), v(1:n,1:n), err = err)
3056 else
3057 call eigen(a(1:n,1:n), vals(1:n))
3058 end if
3059 if (err%has_error_occurred()) then
3060 flag = err%get_error_flag()
3061 return
3062 end if
3063 end function
3064
3065! ------------------------------------------------------------------------------
3083 function la_sort_eigen(ascend, n, vals, vecs, ldv) &
3084 bind(C, name = "la_sort_eigen") result(flag)
3085 ! Arguments
3086 logical(c_bool), intent(in), value :: ascend
3087 integer(c_int), intent(in), value :: n, ldv
3088 real(c_double), intent(inout) :: vals(*), vecs(ldv,*)
3089 integer(c_int) :: flag
3090
3091 ! Local Variables
3092 type(errors) :: err
3093
3094 ! Error Checking
3095 call err%set_exit_on_error(.false.)
3096 flag = la_no_error
3097 if (ldv < n) then
3098 flag = la_invalid_input_error
3099 return
3100 end if
3101
3102 ! Process
3103 call sort(vals(1:n), vecs(1:n,1:n), logical(ascend), err = err)
3104 if (err%has_error_occurred()) then
3105 flag = err%get_error_flag()
3106 return
3107 end if
3108 end function
3109
3110! ------------------------------------------------------------------------------
3128 function la_sort_eigen_cmplx(ascend, n, vals, vecs, ldv) &
3129 bind(C, name = "la_sort_eigen_cmplx") result(flag)
3130 ! Arguments
3131 logical(c_bool), intent(in), value :: ascend
3132 integer(c_int), intent(in), value :: n, ldv
3133 complex(c_double), intent(inout) :: vals(*), vecs(ldv,*)
3134 integer(c_int) :: flag
3135
3136 ! Local Variables
3137 type(errors) :: err
3138
3139 ! Error Checking
3140 call err%set_exit_on_error(.false.)
3141 flag = la_no_error
3142 if (ldv < n) then
3143 flag = la_invalid_input_error
3144 return
3145 end if
3146
3147 ! Process
3148 call sort(vals(1:n), vecs(1:n,1:n), logical(ascend), err = err)
3149 if (err%has_error_occurred()) then
3150 flag = err%get_error_flag()
3151 return
3152 end if
3153 end function
3154
3155! ------------------------------------------------------------------------------
3156
3157! ------------------------------------------------------------------------------
3158
3159! ------------------------------------------------------------------------------
3160end module
Computes the Cholesky factorization of a symmetric, positive definite matrix.
Definition: linalg.f90:1433
Computes the rank 1 downdate to a Cholesky factored matrix (upper triangular).
Definition: linalg.f90:1639
Computes the rank 1 update to a Cholesky factored matrix (upper triangular).
Definition: linalg.f90:1532
Computes the determinant of a square matrix.
Definition: linalg.f90:434
Multiplies a diagonal matrix with another matrix or array.
Definition: linalg.f90:329
Computes the eigenvalues, and optionally the eigenvectors, of a matrix.
Definition: linalg.f90:3098
Extracts the L and U matrices from the condensed [L\U] storage format used by the lu_factor.
Definition: linalg.f90:717
Forms the full M-by-M orthogonal matrix Q from the elementary reflectors returned by the base QR fact...
Definition: linalg.f90:1031
Computes the LU factorization of an M-by-N matrix.
Definition: linalg.f90:595
Computes the inverse of a square matrix.
Definition: linalg.f90:2778
Computes the Moore-Penrose pseudo-inverse of a M-by-N matrix using the singular value decomposition o...
Definition: linalg.f90:2884
Computes the rank of a matrix.
Definition: linalg.f90:401
Multiplies a general matrix by the orthogonal matrix Q from a QR factorization.
Definition: linalg.f90:1185
Computes the QR factorization of an M-by-N matrix.
Definition: linalg.f90:871
Computes the rank 1 update to an M-by-N QR factored matrix A (M >= N) where , and such that .
Definition: linalg.f90:1334
Performs the rank-1 update to matrix A such that: , where is an M-by-N matrix, is a scalar,...
Definition: linalg.f90:194
Solves a system of Cholesky factored equations.
Definition: linalg.f90:2390
Solves the overdetermined or underdetermined system of M equations of N unknowns....
Definition: linalg.f90:2480
Solves a system of LU-factored equations.
Definition: linalg.f90:2149
Solves a system of M QR-factored equations of N unknowns.
Definition: linalg.f90:2284
Solves a triangular system of equations.
Definition: linalg.f90:2061
Sorts an array.
Definition: linalg.f90:3181
Computes the singular value decomposition of a matrix A. The SVD is defined as: , where is an M-by-M...
Definition: linalg.f90:1927
Computes the trace of a matrix (the sum of the main diagonal elements).
Definition: linalg.f90:353
Computes the triangular matrix operation: , or , where A is a triangular matrix.
Definition: linalg.f90:509
Provides a C-friendly API to the LINALG library. Notice, all C-API LINALG routines begin with the pre...
Definition: linalg_c_api.f90:5
Provides a set of common linear algebra routines.
Definition: linalg.f90:15