29 function la_rank1_update(m, n, alpha, x, y, a, lda) &
30 bind(C, name = "la_rank1_update") result(flag)
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
43 flag = la_invalid_input_error
68 function la_rank1_update_cmplx(m, n, alpha, x, y, a, lda) &
69 bind(C, name = "la_rank1_update_cmplx") result(flag)
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
82 flag = la_invalid_input_error
103 function la_trace(m, n, a, lda, rst)
bind(C, name = "la_trace") &
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
116 flag = la_invalid_input_error
121 rst =
trace(a(1:m,1:n))
137 function la_trace_cmplx(m, n, a, lda, rst) &
138 bind(C, name = "la_trace_cmplx") result(flag)
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
150 flag = la_invalid_input_error
155 rst =
trace(a(1:m,1:n))
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)
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
195 integer(c_int) :: nrowa, nrowb
218 if (lda < nrowa .or. ldb < nrowb .or. ldc < m)
then
219 flag = la_invalid_input_error
224 call dgemm(ta, tb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
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)
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
265 integer(c_int) :: nrowa, nrowb
269 if (opa == la_transpose)
then
271 else if (opa == la_hermitian_transpose)
then
277 if (opb == la_transpose)
then
279 else if (opb == la_hermitian_transpose)
then
285 if (opa == la_transpose .or. opa == la_hermitian_transpose)
then
291 if (opb == la_transpose .or. opb == la_hermitian_transpose)
then
298 if (lda < nrowa .or. ldb < nrowb .or. ldc < m)
then
299 flag = la_invalid_input_error
304 call zgemm(ta, tb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
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)
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
350 integer(c_int) :: nrows, ncols, p
355 call err%set_exit_on_error(.false.)
357 if (lside .and. transb)
then
363 else if (lside .and. .not. transb)
then
369 else if (.not. lside .and. transb)
then
384 if (ldb < nrows .or. ldc < m)
then
385 flag = la_invalid_input_error
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()
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)
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
440 integer(c_int) :: nrows, ncols, p
445 call err%set_exit_on_error(.false.)
448 if (opb == la_transpose .or. opb == la_hermitian_transpose) tb = .true.
449 if (lside .and. tb)
then
454 else if (lside .and. .not. tb)
then
459 else if (.not. lside .and. tb)
then
472 if (ldb < nrows .or. ldc < m)
then
473 flag = la_invalid_input_error
478 call diag_mtx_mult(ls, opb, alpha, a(1:p), b(1:nrows,1:ncols), &
480 if (err%has_error_occurred()) flag = err%get_error_flag()
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") &
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
528 integer(c_int) :: nrows, ncols, p
533 call err%set_exit_on_error(.false.)
536 if (opb == la_transpose .or. opb == la_hermitian_transpose) tb = .true.
537 if (lside .and. tb)
then
542 else if (lside .and. .not. tb)
then
547 else if (.not. lside .and. tb)
then
560 if (ldb < nrows .or. ldc < m)
then
561 flag = la_invalid_input_error
566 call diag_mtx_mult(ls, opb, alpha, a(1:p), b(1:nrows,1:ncols), &
568 if (err%has_error_occurred()) flag = err%get_error_flag()
588 function la_rank(m, n, a, lda, rnk)
bind(C, name="la_rank") result(flag)
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
599 call err%set_exit_on_error(.false.)
602 flag = la_invalid_input_error
607 rnk =
mtx_rank(a(1:m,1:n), err =err)
608 if (err%has_error_occurred()) flag = err%get_error_flag()
628 function la_rank_cmplx(m, n, a, lda, rnk)
bind(C, name="la_rank_cmplx") &
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
640 call err%set_exit_on_error(.false.)
643 flag = la_invalid_input_error
648 rnk =
mtx_rank(a(1:m,1:n), err = err)
649 if (err%has_error_occurred()) flag = err%get_error_flag()
666 function la_det(n, a, lda, d)
bind(C, name="la_det") result(flag)
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
677 call err%set_exit_on_error(.false.)
680 flag = la_invalid_input_error
685 d =
det(a(1:n,1:n), err = err)
686 if (err%has_error_occurred()) flag = err%get_error_flag()
703 function la_det_cmplx(n, a, lda, d)
bind(C, name="la_det_cmplx") result(flag)
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
714 call err%set_exit_on_error(.false.)
717 flag = la_invalid_input_error
722 d =
det(a(1:n,1:n), err = err)
723 if (err%has_error_occurred()) flag = err%get_error_flag()
750 function la_tri_mtx_mult(upper, alpha, n, a, lda, beta, b, ldb) &
751 bind(C, name = "la_tri_mtx_mult") result(flag)
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
762 if (lda < n .or. ldb < n)
then
763 flag = la_invalid_input_error
768 call tri_mtx_mult(
logical(upper), alpha, a(1:n,1:n), beta, b(1:n,1:n))
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)
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
807 if (lda < n .or. ldb < n)
then
808 flag = la_invalid_input_error
813 call tri_mtx_mult(
logical(upper), alpha, a(1:n,1:n), beta, b(1:n,1:n))
834 function la_lu_factor(m, n, a, lda, ipvt)
bind(C, name = "la_lu_factor") &
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
847 call err%set_exit_on_error(.false.)
850 flag = la_invalid_input_error
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()
881 function la_lu_factor_cmplx(m, n, a, lda, ipvt) &
882 bind(C, name = "la_lu_factor_cmplx") result(flag)
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
894 call err%set_exit_on_error(.false.)
897 flag = la_invalid_input_error
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()
930 function la_form_lu(n, a, lda, ipvt, u, ldu, p, ldp) &
931 bind(C, name = "la_form_lu") result(flag)
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
941 if (lda < n .or. ldu < n .or. ldp < n)
then
942 flag = la_invalid_input_error
947 call form_lu(a(1:n,1:n), ipvt(1:n), u(1:n,1:n), p(1:n,1:n))
970 function la_form_lu_cmplx(n, a, lda, ipvt, u, ldu, p, ldp) &
971 bind(C, name = "la_form_lu_cmplx") result(flag)
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
982 if (lda < n .or. ldu < n .or. ldp < n)
then
983 flag = la_invalid_input_error
988 call form_lu(a(1:n,1:n), ipvt(1:n), u(1:n,1:n), p(1:n,1:n))
1011 function la_qr_factor(m, n, a, lda, tau)
bind(C, name = "la_qr_factor") &
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
1021 integer(c_int) :: mn
1024 call err%set_exit_on_error(.false.)
1027 flag = la_invalid_input_error
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()
1060 function la_qr_factor_cmplx(m, n, a, lda, tau) &
1061 bind(C, name = "la_qr_factor_cmplx") result(flag)
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
1070 integer(c_int) :: mn
1073 call err%set_exit_on_error(.false.)
1076 flag = la_invalid_input_error
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()
1113 function la_qr_factor_pvt(m, n, a, lda, tau, jpvt) &
1114 bind(C, name = "la_qr_factor_pvt") result(flag)
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
1124 integer(c_int) :: mn
1127 call err%set_exit_on_error(.false.)
1130 flag = la_invalid_input_error
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()
1167 function la_qr_factor_cmplx_pvt(m, n, a, lda, tau, jpvt) &
1168 bind(C, name = "la_qr_factor_cmplx_pvt") result(flag)
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
1178 integer(c_int) :: mn
1181 call err%set_exit_on_error(.false.)
1184 flag = la_invalid_input_error
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()
1222 function la_form_qr(fullq, m, n, r, ldr, tau, q, ldq) &
1223 bind(C, name = "la_form_qr") result(flag)
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
1234 integer(c_int) :: mn, nq
1237 call err%set_exit_on_error(.false.)
1239 if (ldr < m .or. ldq < m)
then
1240 flag = la_invalid_input_error
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()
1280 function la_form_qr_cmplx(fullq, m, n, r, ldr, tau, q, ldq) &
1281 bind(C, name = "la_form_qr_cmplx") result(flag)
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
1292 integer(c_int) :: mn, nq
1295 call err%set_exit_on_error(.false.)
1297 if (ldr < m .or. ldq < m)
then
1298 flag = la_invalid_input_error
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()
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)
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
1357 integer(c_int) :: mn, nq
1360 call err%set_exit_on_error(.false.)
1362 if (ldr < m .or. ldq < m .or. ldp < n)
then
1363 flag = la_invalid_input_error
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), &
1373 if (err%has_error_occurred())
then
1374 flag = err%get_error_flag()
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)
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
1423 integer(c_int) :: mn, nq
1426 call err%set_exit_on_error(.false.)
1428 if (ldr < m .or. ldq < m .or. ldp < n)
then
1429 flag = la_invalid_input_error
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), &
1439 if (err%has_error_occurred())
then
1440 flag = err%get_error_flag()
1473 function la_mult_qr(lside, trans, m, n, k, a, lda, tau, c, ldc) &
1474 bind(C, name = "la_mult_qr") result(flag)
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
1484 integer(c_int) :: ma, na
1496 call err%set_exit_on_error(.false.)
1498 if (lda < ma .or. ldc < m)
then
1499 flag = la_invalid_input_error
1502 if (k > na .or. k < 0)
then
1503 flag = la_invalid_input_error
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()
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)
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
1555 integer(c_int) :: ma, na
1567 call err%set_exit_on_error(.false.)
1569 if (lda < ma .or. ldc < m)
then
1570 flag = la_invalid_input_error
1573 if (k > na .or. k < 0)
then
1574 flag = la_invalid_input_error
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()
1609 function la_qr_rank1_update(m, n, q, ldq, r, ldr, u, v) &
1610 bind(C, name = "la_qr_rank1_update") result(flag)
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
1618 integer(c_int) :: mn
1621 call err%set_exit_on_error(.false.)
1623 if (ldq < m .or. ldr < m)
then
1624 flag = la_invalid_input_error
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()
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)
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
1668 integer(c_int) :: mn
1671 call err%set_exit_on_error(.false.)
1673 if (ldq < m .or. ldr < m)
then
1674 flag = la_invalid_input_error
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()
1704 function la_cholesky_factor(upper, n, a, lda) &
1705 bind(C, name = "la_cholesky_factor") result(flag)
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
1716 call err%set_exit_on_error(.false.)
1719 flag = la_invalid_input_error
1725 if (err%has_error_occurred())
then
1726 flag = err%get_error_flag()
1748 function la_cholesky_factor_cmplx(upper, n, a, lda) &
1749 bind(C, name = "la_cholesky_factor_cmplx") result(flag)
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
1760 call err%set_exit_on_error(.false.)
1763 flag = la_invalid_input_error
1769 if (err%has_error_occurred())
then
1770 flag = err%get_error_flag()
1791 function la_cholesky_rank1_update(n, r, ldr, u) &
1792 bind(C, name = "la_cholesky_rank1_update") result(flag)
1794 integer(c_int),
intent(in),
value :: n, ldr
1795 real(c_double),
intent(inout) :: r(ldr,*), u(*)
1796 integer(c_int) :: flag
1802 call err%set_exit_on_error(.false.)
1805 flag = la_invalid_input_error
1811 if (err%has_error_occurred())
then
1812 flag = err%get_error_flag()
1833 function la_cholesky_rank1_update_cmplx(n, r, ldr, u) &
1834 bind(C, name = "la_cholesky_rank1_update_cmplx") result(flag)
1836 integer(c_int),
intent(in),
value :: n, ldr
1837 complex(c_double),
intent(inout) :: r(ldr,*), u(*)
1838 integer(c_int) :: flag
1844 call err%set_exit_on_error(.false.)
1847 flag = la_invalid_input_error
1853 if (err%has_error_occurred())
then
1854 flag = err%get_error_flag()
1877 function la_cholesky_rank1_downdate(n, r, ldr, u) &
1878 bind(C, name = "la_cholesky_rank1_downdate") result(flag)
1880 integer(c_int),
intent(in),
value :: n, ldr
1881 real(c_double),
intent(inout) :: r(ldr,*), u(*)
1882 integer(c_int) :: flag
1888 call err%set_exit_on_error(.false.)
1891 flag = la_invalid_input_error
1897 if (err%has_error_occurred())
then
1898 flag = err%get_error_flag()
1921 function la_cholesky_rank1_downdate_cmplx(n, r, ldr, u) &
1922 bind(C, name = "la_cholesky_rank1_downdate_cmplx") result(flag)
1924 integer(c_int),
intent(in),
value :: n, ldr
1925 complex(c_double),
intent(inout) :: r(ldr,*), u(*)
1926 integer(c_int) :: flag
1932 call err%set_exit_on_error(.false.)
1935 flag = la_invalid_input_error
1941 if (err%has_error_occurred())
then
1942 flag = err%get_error_flag()
1975 function la_svd(m, n, a, lda, s, u, ldu, vt, ldv) &
1976 bind(C, name = "la_svd") result(flag)
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
1985 integer(c_int) :: mn
1988 call err%set_exit_on_error(.false.)
1990 if (lda < m .or. ldu < m .or. ldv < n)
then
1991 flag = la_invalid_input_error
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()
2032 function la_svd_cmplx(m, n, a, lda, s, u, ldu, vt, ldv) &
2033 bind(C, name = "la_svd_cmplx") result(flag)
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
2043 integer(c_int) :: mn
2046 call err%set_exit_on_error(.false.)
2048 if (lda < m .or. ldu < m .or. ldv < n)
then
2049 flag = la_invalid_input_error
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()
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)
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
2101 integer(c_int) :: ma
2111 call err%set_exit_on_error(.false.)
2113 if (lda < ma .or. ldb < m)
then
2114 flag = la_invalid_input_error
2120 logical(trans),
logical(nounit), alpha, a(1:ma,1:ma), b(1:m,1:n))
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)
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
2163 integer(c_int) :: ma
2173 call err%set_exit_on_error(.false.)
2175 if (lda < ma .or. ldb < m)
then
2176 flag = la_invalid_input_error
2182 logical(trans),
logical(nounit), alpha, a(1:ma,1:ma), b(1:m,1:n))
2200 function la_solve_lu(m, n, a, lda, ipvt, b, ldb) &
2201 bind(C, name = "la_solve_lu") result(flag)
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
2213 call err%set_exit_on_error(.false.)
2215 if (lda < m .or. ldb < m)
then
2216 flag = la_invalid_input_error
2221 call solve_lu(a(1:m,1:m), ipvt(1:m), b(1:m,1:n))
2239 function la_solve_lu_cmplx(m, n, a, lda, ipvt, b, ldb) &
2240 bind(C, name = "la_solve_lu_cmplx") result(flag)
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
2252 call err%set_exit_on_error(.false.)
2254 if (lda < m .or. ldb < m)
then
2255 flag = la_invalid_input_error
2260 call solve_lu(a(1:m,1:m), ipvt(1:m), b(1:m,1:n))
2285 function la_solve_qr(m, n, k, a, lda, tau, b, ldb) &
2286 bind(C, name = "la_solve_qr") result(flag)
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
2295 integer(c_int) :: minmn
2298 call err%set_exit_on_error(.false.)
2300 if (lda < m .or. ldb < m .or. m < n)
then
2301 flag = la_invalid_input_error
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()
2336 function la_solve_qr_cmplx(m, n, k, a, lda, tau, b, ldb) &
2337 bind(C, name = "la_solve_qr_cmplx") result(flag)
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
2346 integer(c_int) :: minmn
2349 call err%set_exit_on_error(.false.)
2351 if (lda < m .or. ldb < m .or. m < n)
then
2352 flag = la_invalid_input_error
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()
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)
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
2398 integer(c_int) :: minmn, maxmn
2403 call err%set_exit_on_error(.false.)
2405 if (lda < m .or. ldb < maxmn)
then
2406 flag = la_invalid_input_error
2411 call solve_qr(a(1:m,1:n), tau(1:minmn), jpvt(1:n), b(1:maxmn,1:k), &
2413 if (err%has_error_occurred())
then
2414 flag = err%get_error_flag()
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)
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
2452 integer(c_int) :: minmn, maxmn
2457 call err%set_exit_on_error(.false.)
2459 if (lda < m .or. ldb < maxmn)
then
2460 flag = la_invalid_input_error
2465 call solve_qr(a(1:m,1:n), tau(1:minmn), jpvt(1:n), b(1:maxmn,1:k), &
2467 if (err%has_error_occurred())
then
2468 flag = err%get_error_flag()
2490 function la_solve_cholesky(upper, m, n, a, lda, b, ldb) &
2491 bind(C, name = "la_solve_cholesky") result(flag)
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
2503 call err%set_exit_on_error(.false.)
2505 if (lda < m .or. ldb < m)
then
2506 flag = la_invalid_input_error
2531 function la_solve_cholesky_cmplx(upper, m, n, a, lda, b, ldb) &
2532 bind(C, name = "la_solve_cholesky_cmplx") result(flag)
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
2544 call err%set_exit_on_error(.false.)
2546 if (lda < m .or. ldb < m)
then
2547 flag = la_invalid_input_error
2579 function la_solve_least_squares(m, n, k, a, lda, b, ldb) &
2580 bind(C, name = "la_solve_least_squares") result(flag)
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
2588 integer(c_int) :: maxmn
2592 call err%set_exit_on_error(.false.)
2594 if (lda < m .or. ldb < maxmn)
then
2595 flag = la_invalid_input_error
2601 if (err%has_error_occurred())
then
2602 flag = err%get_error_flag()
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)
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
2640 integer(c_int) :: maxmn
2644 call err%set_exit_on_error(.false.)
2646 if (lda < m .or. ldb < maxmn)
then
2647 flag = la_invalid_input_error
2653 if (err%has_error_occurred())
then
2654 flag = err%get_error_flag()
2671 function la_inverse(n, a, lda)
bind(C, name = "la_inverse") result(flag)
2673 integer(c_int),
intent(in),
value :: n, lda
2674 real(c_double),
intent(inout) :: a(lda,*)
2675 integer(c_int) :: flag
2681 call err%set_exit_on_error(.false.)
2684 flag = la_invalid_input_error
2690 if (err%has_error_occurred())
then
2691 flag = err%get_error_flag()
2708 function la_inverse_cmplx(n, a, lda)
bind(C, name = "la_inverse_cmplx") &
2711 integer(c_int),
intent(in),
value :: n, lda
2712 complex(c_double),
intent(inout) :: a(lda,*)
2713 integer(c_int) :: flag
2719 call err%set_exit_on_error(.false.)
2722 flag = la_invalid_input_error
2728 if (err%has_error_occurred())
then
2729 flag = err%get_error_flag()
2750 function la_pinverse(m, n, a, lda, ainv, ldai) &
2751 bind(C, name = "la_pinverse") result(flag)
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
2762 call err%set_exit_on_error(.false.)
2764 if (lda < m .or. ldai < n)
then
2765 flag = la_invalid_input_error
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()
2793 function la_pinverse_cmplx(m, n, a, lda, ainv, ldai) &
2794 bind(C, name = "la_pinverse_cmplx") result(flag)
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
2805 call err%set_exit_on_error(.false.)
2807 if (lda < m .or. ldai < n)
then
2808 flag = la_invalid_input_error
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()
2842 function la_eigen_symm(vecs, n, a, lda, vals) &
2843 bind(C, name = "la_eigen_symm") result(flag)
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
2855 call err%set_exit_on_error(.false.)
2858 flag = la_invalid_input_error
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()
2891 function la_eigen_asymm(vecs, n, a, lda, vals, v, ldv) &
2892 bind(C, name = "la_eigen_asymm") result(flag)
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
2904 call err%set_exit_on_error(.false.)
2907 if (lda < n .or. ldv < n)
then
2908 flag = la_invalid_input_error
2913 flag = la_invalid_input_error
2920 call eigen(a(1:n,1:n), vals(1:n), v(1:n,1:n), err = err)
2922 call eigen(a(1:n,1:n), vals(1:n))
2924 if (err%has_error_occurred())
then
2925 flag = err%get_error_flag()
2964 function la_eigen_gen(vecs, n, a, lda, b, ldb, alpha, beta, v, ldv) &
2965 bind(C, name = "la_eigen_gen") result(flag)
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
2978 call err%set_exit_on_error(.false.)
2981 if (lda < n .or. ldb < n .or. ldv < n)
then
2982 flag = la_invalid_input_error
2986 if (lda < n .or. ldb < n)
then
2987 flag = la_invalid_input_error
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)
2997 call eigen(a(1:n,1:n), b(1:n,1:n), alpha(1:n), beta(1:n), err = err)
2999 if (err%has_error_occurred())
then
3000 flag = err%get_error_flag()
3026 function la_eigen_cmplx(vecs, n, a, lda, vals, v, ldv) &
3027 bind(C, name = "la_eigen_cmplx") result(flag)
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
3039 call err%set_exit_on_error(.false.)
3042 if (lda < n .or. ldv < n)
then
3043 flag = la_invalid_input_error
3048 flag = la_invalid_input_error
3055 call eigen(a(1:n,1:n), vals(1:n), v(1:n,1:n), err = err)
3057 call eigen(a(1:n,1:n), vals(1:n))
3059 if (err%has_error_occurred())
then
3060 flag = err%get_error_flag()
3083 function la_sort_eigen(ascend, n, vals, vecs, ldv) &
3084 bind(C, name = "la_sort_eigen") result(flag)
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
3095 call err%set_exit_on_error(.false.)
3098 flag = la_invalid_input_error
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()
3128 function la_sort_eigen_cmplx(ascend, n, vals, vecs, ldv) &
3129 bind(C, name = "la_sort_eigen_cmplx") result(flag)
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
3140 call err%set_exit_on_error(.false.)
3143 flag = la_invalid_input_error
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()
Computes the Cholesky factorization of a symmetric, positive definite matrix.
Computes the rank 1 downdate to a Cholesky factored matrix (upper triangular).
Computes the rank 1 update to a Cholesky factored matrix (upper triangular).
Computes the determinant of a square matrix.
Multiplies a diagonal matrix with another matrix or array.
Computes the eigenvalues, and optionally the eigenvectors, of a matrix.
Computes the LU factorization of an M-by-N matrix.
Computes the inverse of a square matrix.
Computes the Moore-Penrose pseudo-inverse of a M-by-N matrix using the singular value decomposition o...
Computes the rank of a matrix.
Multiplies a general matrix by the orthogonal matrix Q from a QR factorization.
Computes the QR factorization of an M-by-N matrix.
Computes the rank 1 update to an M-by-N QR factored matrix A (M >= N) where , and such that .
Performs the rank-1 update to matrix A such that: , where is an M-by-N matrix, is a scalar,...
Solves a system of Cholesky factored equations.
Solves the overdetermined or underdetermined system of M equations of N unknowns....
Solves a system of LU-factored equations.
Solves a system of M QR-factored equations of N unknowns.
Solves a triangular system of equations.
Computes the singular value decomposition of a matrix A. The SVD is defined as: , where is an M-by-M...
Computes the trace of a matrix (the sum of the main diagonal elements).
Computes the triangular matrix operation: , or , where A is a triangular matrix.
Provides a C-friendly API to the LINALG library. Notice, all C-API LINALG routines begin with the pre...
Provides a set of common linear algebra routines.