7submodule(
linalg) linalg_factor
12 module subroutine lu_factor_dbl(a, ipvt, err)
14 real(real64),
intent(inout),
dimension(:,:) :: a
15 integer(int32),
intent(out),
dimension(:) :: ipvt
16 class(errors),
intent(inout),
optional,
target :: err
19 integer(int32) :: m, n, mn, flag
20 class(errors),
pointer :: errmgr
21 type(errors),
target :: deferr
22 character(len = 128) :: errmsg
28 if (
present(err))
then
36 if (
size(ipvt) /= mn)
then
38 call errmgr%report_error(
"lu_factor_dbl", &
39 "Incorrectly sized input array IPVT, argument 2.", &
45 call dgetrf(m, n, a, m, ipvt, flag)
53 "Singular matrix encountered (row ", flag,
")"
54 call errmgr%report_warning(
"lu_factor_dbl", trim(errmsg), &
55 la_singular_matrix_error)
63 module subroutine lu_factor_cmplx(a, ipvt, err)
65 complex(real64),
intent(inout),
dimension(:,:) :: a
66 integer(int32),
intent(out),
dimension(:) :: ipvt
67 class(errors),
intent(inout),
optional,
target :: err
70 integer(int32) :: m, n, mn, flag
71 class(errors),
pointer :: errmgr
72 type(errors),
target :: deferr
73 character(len = 128) :: errmsg
79 if (
present(err))
then
87 if (
size(ipvt) /= mn)
then
89 call errmgr%report_error(
"lu_factor_cmplx", &
90 "Incorrectly sized input array IPVT, argument 2.", &
96 call zgetrf(m, n, a, m, ipvt, flag)
104 "Singular matrix encountered (row ", flag,
")"
105 call errmgr%report_warning(
"lu_factor_cmplx", trim(errmsg), &
106 la_singular_matrix_error)
114 module subroutine form_lu_all(lu, ipvt, u, p, err)
116 real(real64),
intent(inout),
dimension(:,:) :: lu
117 integer(int32),
intent(in),
dimension(:) :: ipvt
118 real(real64),
intent(out),
dimension(:,:) :: u, p
119 class(errors),
intent(inout),
optional,
target :: err
122 integer(int32) :: j, jp, n, flag
123 class(errors),
pointer :: errmgr
124 type(errors),
target :: deferr
125 character(len = 128) :: errmsg
128 real(real64),
parameter :: zero = 0.0d0
129 real(real64),
parameter :: one = 1.0d0
133 if (
present(err))
then
141 if (
size(lu, 2) /= n)
then
143 else if (
size(ipvt) /= n)
then
145 else if (
size(u, 1) /= n .or.
size(u, 2) /= n)
then
147 else if (
size(p, 1) /= n .or.
size(p, 2) /= n)
then
152 write(errmsg, 100)
"Input number ", flag, &
153 " is not sized correctly."
154 call errmgr%report_error(
"form_lu_all", trim(errmsg), &
160 call dlaset(
'A', n, n, zero, one, p, n)
166 if (j /= jp)
call swap(p(j,1:n), p(jp,1:n))
172 if (j > 1) lu(1:j-1,j) = zero
181 module subroutine form_lu_all_cmplx(lu, ipvt, u, p, err)
183 complex(real64),
intent(inout),
dimension(:,:) :: lu
184 integer(int32),
intent(in),
dimension(:) :: ipvt
185 complex(real64),
intent(out),
dimension(:,:) :: u
186 real(real64),
intent(out),
dimension(:,:) :: p
187 class(errors),
intent(inout),
optional,
target :: err
190 integer(int32) :: j, jp, n, flag
191 class(errors),
pointer :: errmgr
192 type(errors),
target :: deferr
193 character(len = 128) :: errmsg
196 real(real64),
parameter :: zero = 0.0d0
197 real(real64),
parameter :: one = 1.0d0
198 complex(real64),
parameter :: c_zero = (0.0d0, 0.0d0)
199 complex(real64),
parameter :: c_one = (1.0d0, 0.0d0)
203 if (
present(err))
then
211 if (
size(lu, 2) /= n)
then
213 else if (
size(ipvt) /= n)
then
215 else if (
size(u, 1) /= n .or.
size(u, 2) /= n)
then
217 else if (
size(p, 1) /= n .or.
size(p, 2) /= n)
then
222 write(errmsg, 100)
"Input number ", flag, &
223 " is not sized correctly."
224 call errmgr%report_error(
"form_lu_all_cmplx", trim(errmsg), &
230 call dlaset(
'A', n, n, zero, one, p, n)
236 if (j /= jp)
call swap(p(j,1:n), p(jp,1:n))
242 if (j > 1) lu(1:j-1,j) = c_zero
251 module subroutine form_lu_only(lu, u, err)
253 real(real64),
intent(inout),
dimension(:,:) :: lu
254 real(real64),
intent(out),
dimension(:,:) :: u
255 class(errors),
intent(inout),
optional,
target :: err
258 integer(int32) :: j, n, flag
259 class(errors),
pointer :: errmgr
260 type(errors),
target :: deferr
261 character(len = 128) :: errmsg
264 real(real64),
parameter :: zero = 0.0d0
265 real(real64),
parameter :: one = 1.0d0
269 if (
present(err))
then
277 if (
size(lu, 2) /= n)
then
279 else if (
size(u, 1) /= n .or.
size(u, 2) /= n)
then
284 write(errmsg, 100)
"Input number ", flag, &
285 " is not sized correctly."
286 call errmgr%report_error(
"form_lu_only", trim(errmsg), &
297 if (j > 1) lu(1:j-1,j) = zero
306 module subroutine form_lu_only_cmplx(lu, u, err)
308 complex(real64),
intent(inout),
dimension(:,:) :: lu
309 complex(real64),
intent(out),
dimension(:,:) :: u
310 class(errors),
intent(inout),
optional,
target :: err
313 integer(int32) :: j, n, flag
314 class(errors),
pointer :: errmgr
315 type(errors),
target :: deferr
316 character(len = 128) :: errmsg
319 complex(real64),
parameter :: zero = (0.0d0, 0.0d0)
320 complex(real64),
parameter :: one = (1.0d0, 0.0d0)
324 if (
present(err))
then
332 if (
size(lu, 2) /= n)
then
334 else if (
size(u, 1) /= n .or.
size(u, 2) /= n)
then
339 write(errmsg, 100)
"Input number ", flag, &
340 " is not sized correctly."
341 call errmgr%report_error(
"form_lu_only_cmplx", trim(errmsg), &
352 if (j > 1) lu(1:j-1,j) = zero
363 module subroutine qr_factor_no_pivot(a, tau, work, olwork, err)
365 real(real64),
intent(inout),
dimension(:,:) :: a
366 real(real64),
intent(out),
dimension(:) :: tau
367 real(real64),
intent(out),
target,
dimension(:),
optional :: work
368 integer(int32),
intent(out),
optional :: olwork
369 class(errors),
intent(inout),
optional,
target :: err
372 integer(int32) :: m, n, mn, istat, lwork, flag
373 real(real64),
dimension(1) :: temp
374 real(real64),
pointer,
dimension(:) :: wptr
375 real(real64),
allocatable,
target,
dimension(:) :: wrk
376 class(errors),
pointer :: errmgr
377 type(errors),
target :: deferr
383 if (
present(err))
then
390 if (
size(tau) /= mn)
then
392 call errmgr%report_error(
"qr_factor_no_pivot", &
393 "Incorrectly sized input array TAU, argument 2.", &
399 call dgeqrf(m, n, a, m, tau, temp, -1, flag)
400 lwork = int(temp(1), int32)
401 if (
present(olwork))
then
407 if (
present(work))
then
408 if (
size(work) < lwork)
then
410 call errmgr%report_error(
"qr_factor_no_pivot", &
411 "Incorrectly sized input array WORK, argument 3.", &
415 wptr => work(1:lwork)
417 allocate(wrk(lwork), stat = istat)
420 call errmgr%report_error(
"qr_factor_no_pivot", &
421 "Insufficient memory available.", &
422 la_out_of_memory_error)
429 call dgeqrf(m, n, a, m, tau, wptr, lwork, flag)
433 module subroutine qr_factor_no_pivot_cmplx(a, tau, work, olwork, err)
435 complex(real64),
intent(inout),
dimension(:,:) :: a
436 complex(real64),
intent(out),
dimension(:) :: tau
437 complex(real64),
intent(out),
target,
dimension(:),
optional :: work
438 integer(int32),
intent(out),
optional :: olwork
439 class(errors),
intent(inout),
optional,
target :: err
442 integer(int32) :: m, n, mn, istat, lwork, flag
443 complex(real64),
dimension(1) :: temp
444 complex(real64),
pointer,
dimension(:) :: wptr
445 complex(real64),
allocatable,
target,
dimension(:) :: wrk
446 class(errors),
pointer :: errmgr
447 type(errors),
target :: deferr
453 if (
present(err))
then
460 if (
size(tau) /= mn)
then
462 call errmgr%report_error(
"qr_factor_no_pivot_cmplx", &
463 "Incorrectly sized input array TAU, argument 2.", &
469 call zgeqrf(m, n, a, m, tau, temp, -1, flag)
470 lwork = int(temp(1), int32)
471 if (
present(olwork))
then
477 if (
present(work))
then
478 if (
size(work) < lwork)
then
480 call errmgr%report_error(
"qr_factor_no_pivot_cmplx", &
481 "Incorrectly sized input array WORK, argument 3.", &
485 wptr => work(1:lwork)
487 allocate(wrk(lwork), stat = istat)
490 call errmgr%report_error(
"qr_factor_no_pivot_cmplx", &
491 "Insufficient memory available.", &
492 la_out_of_memory_error)
499 call zgeqrf(m, n, a, m, tau, wptr, lwork, flag)
503 module subroutine qr_factor_pivot(a, tau, jpvt, work, olwork, err)
505 real(real64),
intent(inout),
dimension(:,:) :: a
506 real(real64),
intent(out),
dimension(:) :: tau
507 integer(int32),
intent(inout),
dimension(:) :: jpvt
508 real(real64),
intent(out),
target,
dimension(:),
optional :: work
509 integer(int32),
intent(out),
optional :: olwork
510 class(errors),
intent(inout),
optional,
target :: err
513 integer(int32) :: m, n, mn, istat, lwork, flag
514 real(real64),
dimension(1) :: temp
515 real(real64),
pointer,
dimension(:) :: wptr
516 real(real64),
allocatable,
target,
dimension(:) :: wrk
517 class(errors),
pointer :: errmgr
518 type(errors),
target :: deferr
519 character(len = 128) :: errmsg
525 if (
present(err))
then
533 if (
size(tau) /= mn)
then
535 else if (
size(jpvt) /= n)
then
540 write(errmsg, 100)
"Input number ", flag, &
541 " is not sized correctly."
542 call errmgr%report_error(
"qr_factor_pivot", trim(errmsg), &
548 call dgeqp3(m, n, a, m, jpvt, tau, temp, -1, flag)
549 lwork = int(temp(1), int32)
550 if (
present(olwork))
then
556 if (
present(work))
then
557 if (
size(work) < lwork)
then
559 call errmgr%report_error(
"qr_factor_pivot", &
560 "Incorrectly sized input array WORK, argument 4.", &
564 wptr => work(1:lwork)
566 allocate(wrk(lwork), stat = istat)
569 call errmgr%report_error(
"qr_factor_pivot", &
570 "Insufficient memory available.", &
571 la_out_of_memory_error)
578 call dgeqp3(m, n, a, m, jpvt, tau, wptr, lwork, flag)
581 if (
allocated(wrk))
deallocate(wrk)
588 module subroutine qr_factor_pivot_cmplx(a, tau, jpvt, work, olwork, rwork, &
591 complex(real64),
intent(inout),
dimension(:,:) :: a
592 complex(real64),
intent(out),
dimension(:) :: tau
593 integer(int32),
intent(inout),
dimension(:) :: jpvt
594 complex(real64),
intent(out),
target,
dimension(:),
optional :: work
595 integer(int32),
intent(out),
optional :: olwork
596 real(real64),
intent(out),
target,
dimension(:),
optional :: rwork
597 class(errors),
intent(inout),
optional,
target :: err
600 integer(int32) :: m, n, mn, istat, lwork, flag
601 complex(real64),
dimension(1) :: temp
602 complex(real64),
pointer,
dimension(:) :: wptr
603 complex(real64),
allocatable,
target,
dimension(:) :: wrk
604 real(real64),
pointer,
dimension(:) :: rptr
605 real(real64),
allocatable,
target,
dimension(:) :: rwrk
606 class(errors),
pointer :: errmgr
607 type(errors),
target :: deferr
608 character(len = 128) :: errmsg
614 if (
present(err))
then
622 if (
size(tau) /= mn)
then
624 else if (
size(jpvt) /= n)
then
629 write(errmsg, 100)
"Input number ", flag, &
630 " is not sized correctly."
631 call errmgr%report_error(
"qr_factor_pivot_cmplx", trim(errmsg), &
635 if (
present(rwork))
then
636 if (
size(rwork) < 2 * n)
then
637 call errmgr%report_error(
"qr_factor_pivot_cmplx", &
638 "Incorrectly sized input array RWORK, argument 6.", &
644 allocate(rwrk(2 * n), stat = flag)
646 call errmgr%report_error(
"qr_factor_pivot_cmplx", &
647 "Insufficient memory available.", &
648 la_out_of_memory_error)
655 call zgeqp3(m, n, a, m, jpvt, tau, temp, -1, rptr, flag)
656 lwork = int(temp(1), int32)
657 if (
present(olwork))
then
663 if (
present(work))
then
664 if (
size(work) < lwork)
then
666 call errmgr%report_error(
"qr_factor_pivot_cmplx", &
667 "Incorrectly sized input array WORK, argument 4.", &
671 wptr => work(1:lwork)
673 allocate(wrk(lwork), stat = istat)
676 call errmgr%report_error(
"qr_factor_pivot_cmplx", &
677 "Insufficient memory available.", &
678 la_out_of_memory_error)
685 call zgeqp3(m, n, a, m, jpvt, tau, wptr, lwork, rptr, flag)
688 if (
allocated(wrk))
deallocate(wrk)
695 module subroutine form_qr_no_pivot(r, tau, q, work, olwork, err)
697 real(real64),
intent(inout),
dimension(:,:) :: r
698 real(real64),
intent(in),
dimension(:) :: tau
699 real(real64),
intent(out),
dimension(:,:) :: q
700 real(real64),
intent(out),
target,
dimension(:),
optional :: work
701 integer(int32),
intent(out),
optional :: olwork
702 class(errors),
intent(inout),
optional,
target :: err
705 real(real64),
parameter :: zero = 0.0d0
708 integer(int32) :: j, m, n, mn, qcol, istat, flag, lwork
709 real(real64),
pointer,
dimension(:) :: wptr
710 real(real64),
allocatable,
target,
dimension(:) :: wrk
711 real(real64),
dimension(1) :: temp
712 class(errors),
pointer :: errmgr
713 type(errors),
target :: deferr
714 character(len = 128) :: errmsg
721 if (
present(err))
then
729 if (
size(tau) /= mn)
then
731 else if (
size(q, 1) /= m .or. (qcol /= m .and. qcol /= n))
then
733 else if (qcol == n .and. m < n)
then
738 write(errmsg, 100)
"Input number ", flag, &
739 " is not sized correctly."
740 call errmgr%report_error(
"form_qr_no_pivot", trim(errmsg), &
746 call dorgqr(m, qcol, mn, q, m, tau, temp, -1, flag)
747 lwork = int(temp(1), int32)
748 if (
present(olwork))
then
754 if (
present(work))
then
755 if (
size(work) < lwork)
then
757 call errmgr%report_error(
"form_qr_no_pivot", &
758 "Incorrectly sized input array WORK, argument 4.", &
762 wptr => work(1:lwork)
764 allocate(wrk(lwork), stat = istat)
767 call errmgr%report_error(
"form_qr_no_pivot", &
768 "Insufficient memory available.", &
769 la_out_of_memory_error)
778 q(j+1:m,j) = r(j+1:m,j)
783 call dorgqr(m, qcol, mn, q, m, tau, wptr, lwork, flag)
786 if (
allocated(wrk))
deallocate(wrk)
793 module subroutine form_qr_no_pivot_cmplx(r, tau, q, work, olwork, err)
795 complex(real64),
intent(inout),
dimension(:,:) :: r
796 complex(real64),
intent(in),
dimension(:) :: tau
797 complex(real64),
intent(out),
dimension(:,:) :: q
798 complex(real64),
intent(out),
target,
dimension(:),
optional :: work
799 integer(int32),
intent(out),
optional :: olwork
800 class(errors),
intent(inout),
optional,
target :: err
803 complex(real64),
parameter :: zero = (0.0d0, 0.0d0)
806 integer(int32) :: j, m, n, mn, qcol, istat, flag, lwork
807 complex(real64),
pointer,
dimension(:) :: wptr
808 complex(real64),
allocatable,
target,
dimension(:) :: wrk
809 complex(real64),
dimension(1) :: temp
810 class(errors),
pointer :: errmgr
811 type(errors),
target :: deferr
812 character(len = 128) :: errmsg
819 if (
present(err))
then
827 if (
size(tau) /= mn)
then
829 else if (
size(q, 1) /= m .or. (qcol /= m .and. qcol /= n))
then
831 else if (qcol == n .and. m < n)
then
836 write(errmsg, 100)
"Input number ", flag, &
837 " is not sized correctly."
838 call errmgr%report_error(
"form_qr_no_pivot_cmplx", trim(errmsg), &
844 call zungqr(m, qcol, mn, q, m, tau, temp, -1, flag)
845 lwork = int(temp(1), int32)
846 if (
present(olwork))
then
852 if (
present(work))
then
853 if (
size(work) < lwork)
then
855 call errmgr%report_error(
"form_qr_no_pivot_cmplx", &
856 "Incorrectly sized input array WORK, argument 4.", &
860 wptr => work(1:lwork)
862 allocate(wrk(lwork), stat = istat)
865 call errmgr%report_error(
"form_qr_no_pivot_cmplx", &
866 "Insufficient memory available.", &
867 la_out_of_memory_error)
876 q(j+1:m,j) = r(j+1:m,j)
881 call zungqr(m, qcol, mn, q, m, tau, wptr, lwork, flag)
884 if (
allocated(wrk))
deallocate(wrk)
891 module subroutine form_qr_pivot(r, tau, pvt, q, p, work, olwork, err)
893 real(real64),
intent(inout),
dimension(:,:) :: r
894 real(real64),
intent(in),
dimension(:) :: tau
895 integer(int32),
intent(in),
dimension(:) :: pvt
896 real(real64),
intent(out),
dimension(:,:) :: q, p
897 real(real64),
intent(out),
target,
dimension(:),
optional :: work
898 integer(int32),
intent(out),
optional :: olwork
899 class(errors),
intent(inout),
optional,
target :: err
902 real(real64),
parameter :: zero = 0.0d0
903 real(real64),
parameter :: one = 1.0d0
906 integer(int32) :: j, jp, m, n, mn, flag
907 class(errors),
pointer :: errmgr
908 type(errors),
target :: deferr
909 character(len = 128) :: errmsg
915 if (
present(err))
then
923 if (
size(tau) /= mn)
then
925 else if (
size(pvt) /= n)
then
927 else if (
size(q, 1) /= m .or. &
928 (
size(q, 2) /= m .and.
size(q, 2) /= n))
then
930 else if (
size(q, 2) == n .and. m < n)
then
932 else if (
size(p, 1) /= n .or.
size(p, 2) /= n)
then
937 write(errmsg, 100)
"Input number ", flag, &
938 " is not sized correctly."
939 call errmgr%report_error(
"form_qr_pivot", trim(errmsg), &
945 call form_qr_no_pivot(r, tau, q, work = work, olwork = olwork, &
947 if (
present(olwork))
return
948 if (errmgr%has_error_occurred())
return
962 module subroutine form_qr_pivot_cmplx(r, tau, pvt, q, p, work, olwork, err)
964 complex(real64),
intent(inout),
dimension(:,:) :: r
965 complex(real64),
intent(in),
dimension(:) :: tau
966 integer(int32),
intent(in),
dimension(:) :: pvt
967 complex(real64),
intent(out),
dimension(:,:) :: q, p
968 complex(real64),
intent(out),
target,
dimension(:),
optional :: work
969 integer(int32),
intent(out),
optional :: olwork
970 class(errors),
intent(inout),
optional,
target :: err
973 complex(real64),
parameter :: zero = (0.0d0, 0.0d0)
974 complex(real64),
parameter :: one = (1.0d0, 0.0d0)
977 integer(int32) :: j, jp, m, n, mn, flag
978 class(errors),
pointer :: errmgr
979 type(errors),
target :: deferr
980 character(len = 128) :: errmsg
986 if (
present(err))
then
994 if (
size(tau) /= mn)
then
996 else if (
size(pvt) /= n)
then
998 else if (
size(q, 1) /= m .or. &
999 (
size(q, 2) /= m .and.
size(q, 2) /= n))
then
1001 else if (
size(q, 2) == n .and. m < n)
then
1003 else if (
size(p, 1) /= n .or.
size(p, 2) /= n)
then
1008 write(errmsg, 100)
"Input number ", flag, &
1009 " is not sized correctly."
1010 call errmgr%report_error(
"form_qr_pivot_cmplx", trim(errmsg), &
1011 la_array_size_error)
1016 call form_qr_no_pivot_cmplx(r, tau, q, work = work, olwork = olwork, &
1018 if (
present(olwork))
return
1019 if (errmgr%has_error_occurred())
return
1033 module subroutine mult_qr_mtx(lside, trans, a, tau, c, work, olwork, err)
1035 logical,
intent(in) :: lside, trans
1036 real(real64),
intent(in),
dimension(:) :: tau
1037 real(real64),
intent(inout),
dimension(:,:) :: a, c
1038 real(real64),
intent(out),
target,
dimension(:),
optional :: work
1039 integer(int32),
intent(out),
optional :: olwork
1040 class(errors),
intent(inout),
optional,
target :: err
1043 real(real64),
parameter :: one = 1.0d0
1046 character :: side, t
1047 integer(int32) :: m, n, k, nrowa, istat, flag, lwork
1048 real(real64),
pointer,
dimension(:) :: wptr
1049 real(real64),
allocatable,
target,
dimension(:) :: wrk
1050 real(real64),
dimension(1) :: temp
1051 class(errors),
pointer :: errmgr
1052 type(errors),
target :: deferr
1053 character(len = 128) :: errmsg
1071 if (
present(err))
then
1081 if (
size(a, 1) /= m .or.
size(a, 2) < k) flag = 3
1084 if (
size(a, 1) /= n .or.
size(a, 2) < k) flag = 3
1088 write(errmsg, 100)
"Input number ", flag, &
1089 " is not sized correctly."
1090 call errmgr%report_error(
"mult_qr_mtx", trim(errmsg), &
1091 la_array_size_error)
1096 call dormqr(side, t, m, n, k, a, nrowa, tau, c, m, temp, -1, flag)
1097 lwork = int(temp(1), int32)
1098 if (
present(olwork))
then
1104 if (
present(work))
then
1105 if (
size(work) < lwork)
then
1107 call errmgr%report_error(
"mult_qr_mtx", &
1108 "Incorrectly sized input array WORK, argument 6.", &
1109 la_array_size_error)
1112 wptr => work(1:lwork)
1114 allocate(wrk(lwork), stat = istat)
1115 if (istat /= 0)
then
1117 call errmgr%report_error(
"mult_qr_mtx", &
1118 "Insufficient memory available.", &
1119 la_out_of_memory_error)
1126 call dormqr(side, t, m, n, k, a, nrowa, tau, c, m, wptr, lwork, flag)
1133 module subroutine mult_qr_mtx_cmplx(lside, trans, a, tau, c, work, olwork, err)
1135 logical,
intent(in) :: lside, trans
1136 complex(real64),
intent(in),
dimension(:) :: tau
1137 complex(real64),
intent(inout),
dimension(:,:) :: a, c
1138 complex(real64),
intent(out),
target,
dimension(:),
optional :: work
1139 integer(int32),
intent(out),
optional :: olwork
1140 class(errors),
intent(inout),
optional,
target :: err
1143 complex(real64),
parameter :: one = (1.0d0, 0.0d0)
1146 character :: side, t
1147 integer(int32) :: m, n, k, nrowa, istat, flag, lwork
1148 complex(real64),
pointer,
dimension(:) :: wptr
1149 complex(real64),
allocatable,
target,
dimension(:) :: wrk
1150 complex(real64),
dimension(1) :: temp
1151 class(errors),
pointer :: errmgr
1152 type(errors),
target :: deferr
1153 character(len = 128) :: errmsg
1171 if (
present(err))
then
1181 if (
size(a, 1) /= m .or.
size(a, 2) < k) flag = 3
1184 if (
size(a, 1) /= n .or.
size(a, 2) < k) flag = 3
1188 write(errmsg, 100)
"Input number ", flag, &
1189 " is not sized correctly."
1190 call errmgr%report_error(
"mult_qr_mtx_cmplx", trim(errmsg), &
1191 la_array_size_error)
1196 call zunmqr(side, t, m, n, k, a, nrowa, tau, c, m, temp, -1, flag)
1197 lwork = int(temp(1), int32)
1198 if (
present(olwork))
then
1204 if (
present(work))
then
1205 if (
size(work) < lwork)
then
1207 call errmgr%report_error(
"mult_qr_mtx_cmplx", &
1208 "Incorrectly sized input array WORK, argument 6.", &
1209 la_array_size_error)
1212 wptr => work(1:lwork)
1214 allocate(wrk(lwork), stat = istat)
1215 if (istat /= 0)
then
1217 call errmgr%report_error(
"mult_qr_mtx_cmplx", &
1218 "Insufficient memory available.", &
1219 la_out_of_memory_error)
1226 call zunmqr(side, t, m, n, k, a, nrowa, tau, c, m, wptr, lwork, flag)
1233 module subroutine mult_qr_vec(trans, a, tau, c, work, olwork, err)
1235 logical,
intent(in) :: trans
1236 real(real64),
intent(inout),
dimension(:,:) :: a
1237 real(real64),
intent(in),
dimension(:) :: tau
1238 real(real64),
intent(inout),
dimension(:) :: c
1239 real(real64),
intent(out),
target,
dimension(:),
optional :: work
1240 integer(int32),
intent(out),
optional :: olwork
1241 class(errors),
intent(inout),
optional,
target :: err
1244 real(real64),
parameter :: one = 1.0d0
1247 character :: side, t
1248 integer(int32) :: m, k, nrowa, istat, flag, lwork
1249 real(real64),
pointer,
dimension(:) :: wptr
1250 real(real64),
allocatable,
target,
dimension(:) :: wrk
1251 real(real64),
dimension(1) :: temp
1252 class(errors),
pointer :: errmgr
1253 type(errors),
target :: deferr
1254 character(len = 128) :: errmsg
1266 if (
present(err))
then
1274 if (
size(a, 1) /= m .or.
size(a, 2) < k) flag = 3
1277 write(errmsg, 100)
"Input number ", flag, &
1278 " is not sized correctly."
1279 call errmgr%report_error(
"mult_qr_vec", trim(errmsg), &
1280 la_array_size_error)
1285 call dormqr(side, t, m, 1, k, a, nrowa, tau, c, m, temp, -1, flag)
1286 lwork = int(temp(1), int32)
1287 if (
present(olwork))
then
1293 if (
present(work))
then
1294 if (
size(work) < lwork)
then
1296 call errmgr%report_error(
"mult_qr_vec", &
1297 "Incorrectly sized input array WORK, argument 6.", &
1298 la_array_size_error)
1301 wptr => work(1:lwork)
1303 allocate(wrk(lwork), stat = istat)
1304 if (istat /= 0)
then
1306 call errmgr%report_error(
"mult_qr_vec", &
1307 "Insufficient memory available.", &
1308 la_out_of_memory_error)
1315 call dormqr(side, t, m, 1, k, a, nrowa, tau, c, m, wptr, lwork, flag)
1322 module subroutine mult_qr_vec_cmplx(trans, a, tau, c, work, olwork, err)
1324 logical,
intent(in) :: trans
1325 complex(real64),
intent(inout),
dimension(:,:) :: a
1326 complex(real64),
intent(in),
dimension(:) :: tau
1327 complex(real64),
intent(inout),
dimension(:) :: c
1328 complex(real64),
intent(out),
target,
dimension(:),
optional :: work
1329 integer(int32),
intent(out),
optional :: olwork
1330 class(errors),
intent(inout),
optional,
target :: err
1333 complex(real64),
parameter :: one = (1.0d0, 0.0d0)
1336 character :: side, t
1337 integer(int32) :: m, k, nrowa, istat, flag, lwork
1338 complex(real64),
pointer,
dimension(:) :: wptr
1339 complex(real64),
allocatable,
target,
dimension(:) :: wrk
1340 complex(real64),
dimension(1) :: temp
1341 class(errors),
pointer :: errmgr
1342 type(errors),
target :: deferr
1343 character(len = 128) :: errmsg
1355 if (
present(err))
then
1363 if (
size(a, 1) /= m .or.
size(a, 2) < k) flag = 3
1366 write(errmsg, 100)
"Input number ", flag, &
1367 " is not sized correctly."
1368 call errmgr%report_error(
"mult_qr_vec", trim(errmsg), &
1369 la_array_size_error)
1374 call zunmqr(side, t, m, 1, k, a, nrowa, tau, c, m, temp, -1, flag)
1375 lwork = int(temp(1), int32)
1376 if (
present(olwork))
then
1382 if (
present(work))
then
1383 if (
size(work) < lwork)
then
1385 call errmgr%report_error(
"mult_qr_vec", &
1386 "Incorrectly sized input array WORK, argument 6.", &
1387 la_array_size_error)
1390 wptr => work(1:lwork)
1392 allocate(wrk(lwork), stat = istat)
1393 if (istat /= 0)
then
1395 call errmgr%report_error(
"mult_qr_vec", &
1396 "Insufficient memory available.", &
1397 la_out_of_memory_error)
1404 call zunmqr(side, t, m, 1, k, a, nrowa, tau, c, m, wptr, lwork, flag)
1411 module subroutine qr_rank1_update_dbl(q, r, u, v, work, err)
1413 real(real64),
intent(inout),
dimension(:,:) :: q, r
1414 real(real64),
intent(inout),
dimension(:) :: u, v
1415 real(real64),
intent(out),
target,
optional,
dimension(:) :: work
1416 class(errors),
intent(inout),
optional,
target :: err
1420 integer(int32) :: m, n, k, lwork, istat, flag
1421 real(real64),
pointer,
dimension(:) :: wptr
1422 real(real64),
allocatable,
target,
dimension(:) :: wrk
1423 class(errors),
pointer :: errmgr
1424 type(errors),
target :: deferr
1425 character(len = 128) :: errmsg
1431 full =
size(q, 2) == m
1433 if (
present(err))
then
1443 else if (.not.full .and.
size(q, 2) /= k)
then
1445 else if (
size(r, 1) /= m)
then
1447 else if (
size(u) /= m)
then
1449 else if (
size(v) /= n)
then
1454 write(errmsg, 100)
"Input number ", flag, &
1455 " is not sized correctly."
1456 call errmgr%report_error(
"qr_rank1_update", trim(errmsg), &
1457 la_array_size_error)
1462 if (
present(work))
then
1463 if (
size(work) < lwork)
then
1465 call errmgr%report_error(
"qr_rank1_update", &
1466 "Incorrectly sized input array WORK, argument 5.", &
1467 la_array_size_error)
1470 wptr => work(1:lwork)
1472 allocate(wrk(lwork), stat = istat)
1473 if (istat /= 0)
then
1475 call errmgr%report_error(
"qr_rank1_update", &
1476 "Insufficient memory available.", &
1477 la_out_of_memory_error)
1484 call dqr1up(m, n, k, q, m, r, m, u, v, wptr)
1487 if (
allocated(wrk))
deallocate(wrk)
1494 module subroutine qr_rank1_update_cmplx(q, r, u, v, work, rwork, err)
1496 complex(real64),
intent(inout),
dimension(:,:) :: q, r
1497 complex(real64),
intent(inout),
dimension(:) :: u, v
1498 complex(real64),
intent(out),
target,
optional,
dimension(:) :: work
1499 real(real64),
intent(out),
target,
optional,
dimension(:) :: rwork
1500 class(errors),
intent(inout),
optional,
target :: err
1504 integer(int32) :: m, n, k, lwork, istat, flag, lrwork
1505 complex(real64),
pointer,
dimension(:) :: wptr
1506 complex(real64),
allocatable,
target,
dimension(:) :: wrk
1507 real(real64),
pointer,
dimension(:) :: rwptr
1508 real(real64),
allocatable,
target,
dimension(:) :: rwrk
1509 class(errors),
pointer :: errmgr
1510 type(errors),
target :: deferr
1511 character(len = 128) :: errmsg
1517 full =
size(q, 2) == m
1520 if (
present(err))
then
1530 else if (.not.full .and.
size(q, 2) /= k)
then
1532 else if (
size(r, 1) /= m)
then
1534 else if (
size(u) /= m)
then
1536 else if (
size(v) /= n)
then
1541 write(errmsg, 100)
"Input number ", flag, &
1542 " is not sized correctly."
1543 call errmgr%report_error(
"qr_rank1_update_cmplx", trim(errmsg), &
1544 la_array_size_error)
1549 if (
present(work))
then
1550 if (
size(work) < lwork)
then
1552 call errmgr%report_error(
"qr_rank1_update_cmplx", &
1553 "Incorrectly sized input array WORK, argument 5.", &
1554 la_array_size_error)
1557 wptr => work(1:lwork)
1559 allocate(wrk(lwork), stat = istat)
1560 if (istat /= 0)
then
1562 call errmgr%report_error(
"qr_rank1_update_cmplx", &
1563 "Insufficient memory available.", &
1564 la_out_of_memory_error)
1570 if (
present(rwork))
then
1571 if (
size(rwork) < lrwork)
then
1573 call errmgr%report_error(
"qr_rank1_update_cmplx", &
1574 "Incorrectly sized input array RWORK, argument 6.", &
1575 la_array_size_error)
1578 wptr => work(1:lrwork)
1580 allocate(rwrk(lrwork), stat = istat)
1581 if (istat /= 0)
then
1583 call errmgr%report_error(
"qr_rank1_update_cmplx", &
1584 "Insufficient memory available.", &
1585 la_out_of_memory_error)
1592 call zqr1up(m, n, k, q, m, r, m, u, v, wptr, rwptr)
1595 if (
allocated(wrk))
deallocate(wrk)
1604 module subroutine cholesky_factor_dbl(a, upper, err)
1606 real(real64),
intent(inout),
dimension(:,:) :: a
1607 logical,
intent(in),
optional :: upper
1608 class(errors),
intent(inout),
optional,
target :: err
1611 real(real64),
parameter :: zero = 0.0d0
1615 integer(int32) :: i, n, flag
1616 class(errors),
pointer :: errmgr
1617 type(errors),
target :: deferr
1618 character(len = 128) :: errmsg
1622 if (
present(upper))
then
1631 if (
present(err))
then
1638 if (
size(a, 2) /= n)
then
1640 call errmgr%report_error(
"cholesky_factor", &
1641 "The input matrix must be square.", la_array_size_error)
1646 call dpotrf(uplo, n, a, n, flag)
1649 write(errmsg, 100)
"The leading minor of order ", flag, &
1650 " is not positive definite."
1651 call errmgr%report_error(
"cholesky_factor", trim(errmsg), &
1652 la_matrix_format_error)
1656 if (uplo ==
'U')
then
1673 module subroutine cholesky_factor_cmplx(a, upper, err)
1675 complex(real64),
intent(inout),
dimension(:,:) :: a
1676 logical,
intent(in),
optional :: upper
1677 class(errors),
intent(inout),
optional,
target :: err
1680 complex(real64),
parameter :: zero = (0.0d0, 0.0d0)
1684 integer(int32) :: i, n, flag
1685 class(errors),
pointer :: errmgr
1686 type(errors),
target :: deferr
1687 character(len = 128) :: errmsg
1691 if (
present(upper))
then
1700 if (
present(err))
then
1707 if (
size(a, 2) /= n)
then
1709 call errmgr%report_error(
"cholesky_factor_cmplx", &
1710 "The input matrix must be square.", la_array_size_error)
1715 call zpotrf(uplo, n, a, n, flag)
1718 write(errmsg, 100)
"The leading minor of order ", flag, &
1719 " is not positive definite."
1720 call errmgr%report_error(
"cholesky_factor_cmplx", trim(errmsg), &
1721 la_matrix_format_error)
1725 if (uplo ==
'U')
then
1742 module subroutine cholesky_rank1_update_dbl(r, u, work, err)
1744 real(real64),
intent(inout),
dimension(:,:) :: r
1745 real(real64),
intent(inout),
dimension(:) :: u
1746 real(real64),
intent(out),
target,
optional,
dimension(:) :: work
1747 class(errors),
intent(inout),
optional,
target :: err
1750 integer(int32) :: n, lwork, istat, flag
1751 real(real64),
pointer,
dimension(:) :: wptr
1752 real(real64),
allocatable,
target,
dimension(:) :: wrk
1753 class(errors),
pointer :: errmgr
1754 type(errors),
target :: deferr
1755 character(len = 128) :: errmsg
1760 if (
present(err))
then
1768 if (
size(r, 2) /= n)
then
1770 else if (
size(u) /= n)
then
1774 write(errmsg, 100)
"Input number ", flag, &
1775 " is not sized correctly."
1776 call errmgr%report_error(
"cholesky_rank1_update", trim(errmsg), &
1777 la_array_size_error)
1782 if (
present(work))
then
1783 if (
size(work) < lwork)
then
1785 call errmgr%report_error(
"cholesky_rank1_update", &
1786 "The workspace array is too short.", &
1787 la_array_size_error)
1790 wptr => work(1:lwork)
1792 allocate(wrk(lwork), stat = istat)
1793 if (istat /= 0)
then
1794 call errmgr%report_error(
"cholesky_rank1_update", &
1795 "Insufficient memory available.", &
1796 la_out_of_memory_error)
1803 call dch1up(n, r, n, u, wptr)
1810 module subroutine cholesky_rank1_update_cmplx(r, u, work, err)
1812 complex(real64),
intent(inout),
dimension(:,:) :: r
1813 complex(real64),
intent(inout),
dimension(:) :: u
1814 real(real64),
intent(out),
target,
optional,
dimension(:) :: work
1815 class(errors),
intent(inout),
optional,
target :: err
1818 integer(int32) :: n, lwork, istat, flag
1819 real(real64),
pointer,
dimension(:) :: wptr
1820 real(real64),
allocatable,
target,
dimension(:) :: wrk
1821 class(errors),
pointer :: errmgr
1822 type(errors),
target :: deferr
1823 character(len = 128) :: errmsg
1828 if (
present(err))
then
1836 if (
size(r, 2) /= n)
then
1838 else if (
size(u) /= n)
then
1842 write(errmsg, 100)
"Input number ", flag, &
1843 " is not sized correctly."
1844 call errmgr%report_error(
"cholesky_rank1_update_cmplx", &
1846 la_array_size_error)
1851 if (
present(work))
then
1852 if (
size(work) < lwork)
then
1854 call errmgr%report_error(
"cholesky_rank1_update_cmplx", &
1855 "The workspace array is too short.", &
1856 la_array_size_error)
1859 wptr => work(1:lwork)
1861 allocate(wrk(lwork), stat = istat)
1862 if (istat /= 0)
then
1863 call errmgr%report_error(
"cholesky_rank1_update", &
1864 "Insufficient memory available.", &
1865 la_out_of_memory_error)
1872 call zch1up(n, r, n, u, wptr)
1879 module subroutine cholesky_rank1_downdate_dbl(r, u, work, err)
1881 real(real64),
intent(inout),
dimension(:,:) :: r
1882 real(real64),
intent(inout),
dimension(:) :: u
1883 real(real64),
intent(out),
target,
optional,
dimension(:) :: work
1884 class(errors),
intent(inout),
optional,
target :: err
1887 integer(int32) :: n, lwork, istat, flag
1888 real(real64),
pointer,
dimension(:) :: wptr
1889 real(real64),
allocatable,
target,
dimension(:) :: wrk
1890 class(errors),
pointer :: errmgr
1891 type(errors),
target :: deferr
1892 character(len = 128) :: errmsg
1897 if (
present(err))
then
1905 if (
size(r, 2) /= n)
then
1907 else if (
size(u) /= n)
then
1911 write(errmsg, 100)
"Input number ", flag, &
1912 " is not sized correctly."
1913 call errmgr%report_error(
"cholesky_rank1_downdate", trim(errmsg), &
1914 la_array_size_error)
1919 if (
present(work))
then
1920 if (
size(work) < lwork)
then
1922 call errmgr%report_error(
"cholesky_rank1_downdate", &
1923 "The workspace array is too short.", &
1924 la_array_size_error)
1927 wptr => work(1:lwork)
1929 allocate(wrk(lwork), stat = istat)
1930 if (istat /= 0)
then
1931 call errmgr%report_error(
"cholesky_rank1_downdate", &
1932 "Insufficient memory available.", &
1933 la_out_of_memory_error)
1940 call dch1dn(n, r, n, u, wptr, flag)
1943 call errmgr%report_error(
"cholesky_rank1_downdate", &
1944 "The downdated matrix is not positive definite.", &
1945 la_matrix_format_error)
1946 else if (flag == 2)
then
1948 call errmgr%report_error(
"cholesky_rank1_downdate", &
1949 "The input matrix is singular.", la_singular_matrix_error)
1957 module subroutine cholesky_rank1_downdate_cmplx(r, u, work, err)
1959 complex(real64),
intent(inout),
dimension(:,:) :: r
1960 complex(real64),
intent(inout),
dimension(:) :: u
1961 real(real64),
intent(out),
target,
optional,
dimension(:) :: work
1962 class(errors),
intent(inout),
optional,
target :: err
1965 integer(int32) :: n, lwork, istat, flag
1966 real(real64),
pointer,
dimension(:) :: wptr
1967 real(real64),
allocatable,
target,
dimension(:) :: wrk
1968 class(errors),
pointer :: errmgr
1969 type(errors),
target :: deferr
1970 character(len = 128) :: errmsg
1975 if (
present(err))
then
1983 if (
size(r, 2) /= n)
then
1985 else if (
size(u) /= n)
then
1989 write(errmsg, 100)
"Input number ", flag, &
1990 " is not sized correctly."
1991 call errmgr%report_error(
"cholesky_rank1_downdate_cmplx", &
1993 la_array_size_error)
1998 if (
present(work))
then
1999 if (
size(work) < lwork)
then
2001 call errmgr%report_error(
"cholesky_rank1_downdate_cmplx", &
2002 "The workspace array is too short.", &
2003 la_array_size_error)
2006 wptr => work(1:lwork)
2008 allocate(wrk(lwork), stat = istat)
2009 if (istat /= 0)
then
2010 call errmgr%report_error(
"cholesky_rank1_downdate_cmplx", &
2011 "Insufficient memory available.", &
2012 la_out_of_memory_error)
2019 call zch1dn(n, r, n, u, wptr, flag)
2022 call errmgr%report_error(
"cholesky_rank1_downdate_cmplx", &
2023 "The downdated matrix is not positive definite.", &
2024 la_matrix_format_error)
2025 else if (flag == 2)
then
2027 call errmgr%report_error(
"cholesky_rank1_downdate_cmplx", &
2028 "The input matrix is singular.", la_singular_matrix_error)
2038 module subroutine rz_factor_dbl(a, tau, work, olwork, err)
2040 real(real64),
intent(inout),
dimension(:,:) :: a
2041 real(real64),
intent(out),
dimension(:) :: tau
2042 real(real64),
intent(out),
target,
optional,
dimension(:) :: work
2043 integer(int32),
intent(out),
optional :: olwork
2044 class(errors),
intent(inout),
optional,
target :: err
2047 integer(int32) :: m, n, lwork, flag, istat
2048 real(real64),
pointer,
dimension(:) :: wptr
2049 real(real64),
allocatable,
target,
dimension(:) :: wrk
2050 real(real64),
dimension(1) :: temp
2051 class(errors),
pointer :: errmgr
2052 type(errors),
target :: deferr
2053 character(len = 128) :: errmsg
2058 if (
present(err))
then
2066 if (
size(tau) /= m)
then
2071 write(errmsg, 100)
"Input number ", flag, &
2072 " is not sized correctly."
2073 call errmgr%report_error(
"rz_factor", trim(errmsg), &
2074 la_array_size_error)
2079 call dtzrzf(m, n, a, m, tau, temp, -1, flag)
2080 lwork = int(temp(1), int32)
2081 if (
present(olwork))
then
2087 if (
present(work))
then
2088 if (
size(work) < lwork)
then
2090 call errmgr%report_error(
"rz_factor", &
2091 "Incorrectly sized input array WORK, argument 3.", &
2092 la_array_size_error)
2095 wptr => work(1:lwork)
2097 allocate(wrk(lwork), stat = istat)
2098 if (istat /= 0)
then
2100 call errmgr%report_error(
"rz_factor", &
2101 "Insufficient memory available.", &
2102 la_out_of_memory_error)
2109 call dtzrzf(m, n, a, m, tau, wptr, lwork, flag)
2116 module subroutine rz_factor_cmplx(a, tau, work, olwork, err)
2118 complex(real64),
intent(inout),
dimension(:,:) :: a
2119 complex(real64),
intent(out),
dimension(:) :: tau
2120 complex(real64),
intent(out),
target,
optional,
dimension(:) :: work
2121 integer(int32),
intent(out),
optional :: olwork
2122 class(errors),
intent(inout),
optional,
target :: err
2125 integer(int32) :: m, n, lwork, flag, istat
2126 complex(real64),
pointer,
dimension(:) :: wptr
2127 complex(real64),
allocatable,
target,
dimension(:) :: wrk
2128 complex(real64),
dimension(1) :: temp
2129 class(errors),
pointer :: errmgr
2130 type(errors),
target :: deferr
2131 character(len = 128) :: errmsg
2136 if (
present(err))
then
2144 if (
size(tau) /= m)
then
2149 write(errmsg, 100)
"Input number ", flag, &
2150 " is not sized correctly."
2151 call errmgr%report_error(
"rz_factor_cmplx", trim(errmsg), &
2152 la_array_size_error)
2157 call ztzrzf(m, n, a, m, tau, temp, -1, flag)
2158 lwork = int(temp(1), int32)
2159 if (
present(olwork))
then
2165 if (
present(work))
then
2166 if (
size(work) < lwork)
then
2168 call errmgr%report_error(
"rz_factor_cmplx", &
2169 "Incorrectly sized input array WORK, argument 3.", &
2170 la_array_size_error)
2173 wptr => work(1:lwork)
2175 allocate(wrk(lwork), stat = istat)
2176 if (istat /= 0)
then
2178 call errmgr%report_error(
"rz_factor_cmplx", &
2179 "Insufficient memory available.", &
2180 la_out_of_memory_error)
2187 call ztzrzf(m, n, a, m, tau, wptr, lwork, flag)
2194 module subroutine mult_rz_mtx(lside, trans, l, a, tau, c, work, olwork, err)
2196 logical,
intent(in) :: lside, trans
2197 integer(int32),
intent(in) :: l
2198 real(real64),
intent(inout),
dimension(:,:) :: a, c
2199 real(real64),
intent(in),
dimension(:) :: tau
2200 real(real64),
intent(out),
target,
optional,
dimension(:) :: work
2201 integer(int32),
intent(out),
optional :: olwork
2202 class(errors),
intent(inout),
optional,
target :: err
2205 character :: side, t
2206 integer(int32) :: m, n, k, lwork, flag, istat, lda
2207 real(real64),
pointer,
dimension(:) :: wptr
2208 real(real64),
allocatable,
target,
dimension(:) :: wrk
2209 real(real64),
dimension(1) :: temp
2210 class(errors),
pointer :: errmgr
2211 type(errors),
target :: deferr
2212 character(len = 128) :: errmsg
2229 if (
present(err))
then
2238 if (l > m .or. l < 0)
then
2240 else if (k > m)
then
2242 else if (
size(a, 1) < k .or.
size(a, 2) /= m)
then
2246 if (l > n .or. l < 0)
then
2248 else if (k > n)
then
2250 else if (
size(a, 1) < k .or.
size(a, 2) /= n)
then
2256 write(errmsg, 100)
"Input number ", flag, &
2257 " is not sized correctly."
2258 call errmgr%report_error(
"mult_rz_mtx", trim(errmsg), &
2259 la_array_size_error)
2264 call dormrz(side, t, m, n, k, l, a, lda, tau, c, m, temp, -1, flag)
2265 lwork = int(temp(1), int32)
2266 if (
present(olwork))
then
2272 if (
present(work))
then
2273 if (
size(work) < lwork)
then
2275 call errmgr%report_error(
"mult_rz_mtx", &
2276 "Incorrectly sized input array WORK, argument 7.", &
2277 la_array_size_error)
2280 wptr => work(1:lwork)
2282 allocate(wrk(lwork), stat = istat)
2283 if (istat /= 0)
then
2285 call errmgr%report_error(
"mult_rz_mtx", &
2286 "Insufficient memory available.", &
2287 la_out_of_memory_error)
2294 call dormrz(side, t, m, n, k, l, a, lda, tau, c, m, wptr, lwork, flag)
2301 module subroutine mult_rz_mtx_cmplx(lside, trans, l, a, tau, c, work, olwork, err)
2303 logical,
intent(in) :: lside, trans
2304 integer(int32),
intent(in) :: l
2305 complex(real64),
intent(inout),
dimension(:,:) :: a, c
2306 complex(real64),
intent(in),
dimension(:) :: tau
2307 complex(real64),
intent(out),
target,
optional,
dimension(:) :: work
2308 integer(int32),
intent(out),
optional :: olwork
2309 class(errors),
intent(inout),
optional,
target :: err
2312 character :: side, t
2313 integer(int32) :: m, n, k, lwork, flag, istat, lda
2314 complex(real64),
pointer,
dimension(:) :: wptr
2315 complex(real64),
allocatable,
target,
dimension(:) :: wrk
2316 complex(real64),
dimension(1) :: temp
2317 class(errors),
pointer :: errmgr
2318 type(errors),
target :: deferr
2319 character(len = 128) :: errmsg
2336 if (
present(err))
then
2345 if (l > m .or. l < 0)
then
2347 else if (k > m)
then
2349 else if (
size(a, 1) < k .or.
size(a, 2) /= m)
then
2353 if (l > n .or. l < 0)
then
2355 else if (k > n)
then
2357 else if (
size(a, 1) < k .or.
size(a, 2) /= n)
then
2363 write(errmsg, 100)
"Input number ", flag, &
2364 " is not sized correctly."
2365 call errmgr%report_error(
"mult_rz_mtx_cmplx", trim(errmsg), &
2366 la_array_size_error)
2371 call zunmrz(side, t, m, n, k, l, a, lda, tau, c, m, temp, -1, flag)
2372 lwork = int(temp(1), int32)
2373 if (
present(olwork))
then
2379 if (
present(work))
then
2380 if (
size(work) < lwork)
then
2382 call errmgr%report_error(
"mult_rz_mtx_cmplx", &
2383 "Incorrectly sized input array WORK, argument 7.", &
2384 la_array_size_error)
2387 wptr => work(1:lwork)
2389 allocate(wrk(lwork), stat = istat)
2390 if (istat /= 0)
then
2392 call errmgr%report_error(
"mult_rz_mtx_cmplx", &
2393 "Insufficient memory available.", &
2394 la_out_of_memory_error)
2401 call zunmrz(side, t, m, n, k, l, a, lda, tau, c, m, wptr, lwork, flag)
2408 module subroutine mult_rz_vec(trans, l, a, tau, c, work, olwork, err)
2410 logical,
intent(in) :: trans
2411 integer(int32),
intent(in) :: l
2412 real(real64),
intent(inout),
dimension(:,:) :: a
2413 real(real64),
intent(in),
dimension(:) :: tau
2414 real(real64),
intent(inout),
dimension(:) :: c
2415 real(real64),
intent(out),
target,
optional,
dimension(:) :: work
2416 integer(int32),
intent(out),
optional :: olwork
2417 class(errors),
intent(inout),
optional,
target :: err
2420 character :: side, t
2421 integer(int32) :: m, k, lwork, flag, istat, lda
2422 real(real64),
pointer,
dimension(:) :: wptr
2423 real(real64),
allocatable,
target,
dimension(:) :: wrk
2424 real(real64),
dimension(1) :: temp
2425 class(errors),
pointer :: errmgr
2426 type(errors),
target :: deferr
2427 character(len = 128) :: errmsg
2439 if (
present(err))
then
2447 if (l > m .or. l < 0)
then
2449 else if (k > m)
then
2451 else if (
size(a, 1) < k .or.
size(a, 2) /= m)
then
2456 write(errmsg, 100)
"Input number ", flag, &
2457 " is not sized correctly."
2458 call errmgr%report_error(
"mult_rz_vec", trim(errmsg), &
2459 la_array_size_error)
2464 call dormrz(side, t, m, 1, k, l, a, lda, tau, c, m, temp, -1, flag)
2465 lwork = int(temp(1), int32)
2466 if (
present(olwork))
then
2472 if (
present(work))
then
2473 if (
size(work) < lwork)
then
2475 call errmgr%report_error(
"mult_rz_vec", &
2476 "Incorrectly sized input array WORK, argument 6.", &
2477 la_array_size_error)
2480 wptr => work(1:lwork)
2482 allocate(wrk(lwork), stat = istat)
2483 if (istat /= 0)
then
2485 call errmgr%report_error(
"mult_rz_vec", &
2486 "Insufficient memory available.", &
2487 la_out_of_memory_error)
2494 call dormrz(side, t, m, 1, k, l, a, lda, tau, c, m, wptr, lwork, flag)
2501 module subroutine mult_rz_vec_cmplx(trans, l, a, tau, c, work, olwork, err)
2503 logical,
intent(in) :: trans
2504 integer(int32),
intent(in) :: l
2505 complex(real64),
intent(inout),
dimension(:,:) :: a
2506 complex(real64),
intent(in),
dimension(:) :: tau
2507 complex(real64),
intent(inout),
dimension(:) :: c
2508 complex(real64),
intent(out),
target,
optional,
dimension(:) :: work
2509 integer(int32),
intent(out),
optional :: olwork
2510 class(errors),
intent(inout),
optional,
target :: err
2513 character :: side, t
2514 integer(int32) :: m, k, lwork, flag, istat, lda
2515 complex(real64),
pointer,
dimension(:) :: wptr
2516 complex(real64),
allocatable,
target,
dimension(:) :: wrk
2517 complex(real64),
dimension(1) :: temp
2518 class(errors),
pointer :: errmgr
2519 type(errors),
target :: deferr
2520 character(len = 128) :: errmsg
2532 if (
present(err))
then
2540 if (l > m .or. l < 0)
then
2542 else if (k > m)
then
2544 else if (
size(a, 1) < k .or.
size(a, 2) /= m)
then
2549 write(errmsg, 100)
"Input number ", flag, &
2550 " is not sized correctly."
2551 call errmgr%report_error(
"mult_rz_vec_cmplx", trim(errmsg), &
2552 la_array_size_error)
2557 call zunmrz(side, t, m, 1, k, l, a, lda, tau, c, m, temp, -1, flag)
2558 lwork = int(temp(1), int32)
2559 if (
present(olwork))
then
2565 if (
present(work))
then
2566 if (
size(work) < lwork)
then
2568 call errmgr%report_error(
"mult_rz_vec_cmplx", &
2569 "Incorrectly sized input array WORK, argument 6.", &
2570 la_array_size_error)
2573 wptr => work(1:lwork)
2575 allocate(wrk(lwork), stat = istat)
2576 if (istat /= 0)
then
2578 call errmgr%report_error(
"mult_rz_vec_cmplx", &
2579 "Insufficient memory available.", &
2580 la_out_of_memory_error)
2587 call zunmrz(side, t, m, 1, k, l, a, lda, tau, c, m, wptr, lwork, flag)
2596 module subroutine svd_dbl(a, s, u, vt, work, olwork, err)
2598 real(real64),
intent(inout),
dimension(:,:) :: a
2599 real(real64),
intent(out),
dimension(:) :: s
2600 real(real64),
intent(out),
optional,
dimension(:,:) :: u, vt
2601 real(real64),
intent(out),
target,
optional,
dimension(:) :: work
2602 integer(int32),
intent(out),
optional :: olwork
2603 class(errors),
intent(inout),
optional,
target :: err
2606 character :: jobu, jobvt
2607 integer(int32) :: m, n, mn, istat, lwork, flag
2608 real(real64),
pointer,
dimension(:) :: wptr
2609 real(real64),
allocatable,
target,
dimension(:) :: wrk
2610 real(real64),
dimension(1) :: temp
2611 class(errors),
pointer :: errmgr
2612 type(errors),
target :: deferr
2613 character(len = 128) :: errmsg
2619 if (
present(u))
then
2620 if (
size(u, 2) == m)
then
2622 else if (
size(u, 2) == mn)
then
2628 if (
present(vt))
then
2633 if (
present(err))
then
2641 if (
size(s) /= mn)
then
2643 else if (
present(u))
then
2644 if (
size(u, 1) /= m) flag = 3
2645 if (
size(u, 2) /= m .and.
size(u, 2) /= mn) flag = 3
2646 else if (
present(vt))
then
2647 if (
size(vt, 1) /= n .or.
size(vt, 2) /= n) flag = 4
2651 write(errmsg, 100)
"Input number ", flag, &
2652 " is not sized correctly."
2653 call errmgr%report_error(
"svd", trim(errmsg), &
2654 la_array_size_error)
2659 call dgesvd(jobu, jobvt, m, n, a, m, s, temp, m, temp, n, temp, -1, &
2661 lwork = int(temp(1), int32)
2662 if (
present(olwork))
then
2668 if (
present(work))
then
2669 if (
size(work) < lwork)
then
2671 call errmgr%report_error(
"svd", &
2672 "Incorrectly sized input array WORK, argument 5.", &
2673 la_array_size_error)
2676 wptr => work(1:lwork)
2678 allocate(wrk(lwork), stat = istat)
2679 if (istat /= 0)
then
2681 call errmgr%report_error(
"svd", &
2682 "Insufficient memory available.", &
2683 la_out_of_memory_error)
2690 if (
present(u) .and.
present(vt))
then
2691 call dgesvd(jobu, jobvt, m, n, a, m, s, u, m, vt, n, wptr, lwork, &
2693 else if (
present(u) .and. .not.
present(vt))
then
2694 call dgesvd(jobu, jobvt, m, n, a, m, s, u, m, temp, n, wptr, &
2696 else if (.not.
present(u) .and.
present(vt))
then
2697 call dgesvd(jobu, jobvt, m, n, a, m, s, temp, m, vt, n, wptr, &
2700 call dgesvd(jobu, jobvt, m, n, a, m, s, temp, m, temp, n, wptr, &
2706 write(errmsg, 101) flag,
" superdiagonals could not " // &
2707 "converge to zero as part of the QR iteration process."
2708 call errmgr%report_warning(
"svd", errmsg, la_convergence_error)
2717 module subroutine svd_cmplx(a, s, u, vt, work, olwork, rwork, err)
2719 complex(real64),
intent(inout),
dimension(:,:) :: a
2720 real(real64),
intent(out),
dimension(:) :: s
2721 complex(real64),
intent(out),
optional,
dimension(:,:) :: u, vt
2722 complex(real64),
intent(out),
target,
optional,
dimension(:) :: work
2723 integer(int32),
intent(out),
optional :: olwork
2724 real(real64),
intent(out),
target,
optional,
dimension(:) :: rwork
2725 class(errors),
intent(inout),
optional,
target :: err
2728 character :: jobu, jobvt
2729 integer(int32) :: m, n, mn, istat, lwork, flag, lrwork
2730 complex(real64),
pointer,
dimension(:) :: wptr
2731 complex(real64),
allocatable,
target,
dimension(:) :: wrk
2732 complex(real64),
dimension(1) :: temp
2733 real(real64),
dimension(1) :: rtemp
2734 real(real64),
pointer,
dimension(:) :: rwptr
2735 real(real64),
allocatable,
target,
dimension(:) :: rwrk
2736 class(errors),
pointer :: errmgr
2737 type(errors),
target :: deferr
2738 character(len = 128) :: errmsg
2745 if (
present(u))
then
2746 if (
size(u, 2) == m)
then
2748 else if (
size(u, 2) == mn)
then
2754 if (
present(vt))
then
2759 if (
present(err))
then
2767 if (
size(s) /= mn)
then
2769 else if (
present(u))
then
2770 if (
size(u, 1) /= m) flag = 3
2771 if (
size(u, 2) /= m .and.
size(u, 2) /= mn) flag = 3
2772 else if (
present(vt))
then
2773 if (
size(vt, 1) /= n .or.
size(vt, 2) /= n) flag = 4
2777 write(errmsg, 100)
"Input number ", flag, &
2778 " is not sized correctly."
2779 call errmgr%report_error(
"svd_cmplx", trim(errmsg), &
2780 la_array_size_error)
2785 call zgesvd(jobu, jobvt, m, n, a, m, s, temp, m, temp, n, temp, -1, &
2787 lwork = int(temp(1), int32)
2788 if (
present(olwork))
then
2794 if (
present(work))
then
2795 if (
size(work) < lwork)
then
2797 call errmgr%report_error(
"svd_cmplx", &
2798 "Incorrectly sized input array WORK, argument 5.", &
2799 la_array_size_error)
2802 wptr => work(1:lwork)
2804 allocate(wrk(lwork), stat = istat)
2805 if (istat /= 0)
then
2807 call errmgr%report_error(
"svd_cmplx", &
2808 "Insufficient memory available.", &
2809 la_out_of_memory_error)
2815 if (
present(rwork))
then
2816 if (
size(rwork) < lrwork)
then
2818 call errmgr%report_error(
"svd_cmplx", &
2819 "Incorrectly sized input array RWORK, argument 7.", &
2820 la_array_size_error)
2822 rwptr => rwork(1:lrwork)
2824 allocate(rwrk(lrwork), stat = istat)
2825 if (istat /= 0)
then
2827 call errmgr%report_error(
"svd_cmplx", &
2828 "Insufficient memory available.", &
2829 la_out_of_memory_error)
2836 if (
present(u) .and.
present(vt))
then
2837 call zgesvd(jobu, jobvt, m, n, a, m, s, u, m, vt, n, wptr, lwork, &
2839 else if (
present(u) .and. .not.
present(vt))
then
2840 call zgesvd(jobu, jobvt, m, n, a, m, s, u, m, temp, n, wptr, &
2842 else if (.not.
present(u) .and.
present(vt))
then
2843 call zgesvd(jobu, jobvt, m, n, a, m, s, temp, m, vt, n, wptr, &
2846 call zgesvd(jobu, jobvt, m, n, a, m, s, temp, m, temp, n, wptr, &
2852 write(errmsg, 101) flag,
" superdiagonals could not " // &
2853 "converge to zero as part of the QR iteration process."
2854 call errmgr%report_warning(
"svd_cmplx", errmsg, &
2855 la_convergence_error)
2866 module subroutine lq_factor_no_pivot(a, tau, work, olwork, err)
2868 real(real64),
intent(inout),
dimension(:,:) :: a
2869 real(real64),
intent(out),
dimension(:) :: tau
2870 real(real64),
intent(out),
target,
dimension(:),
optional :: work
2871 integer(int32),
intent(out),
optional :: olwork
2872 class(errors),
intent(inout),
optional,
target :: err
2875 integer(int32) :: m, n, mn, istat, lwork, flag
2876 real(real64),
dimension(1) :: temp
2877 real(real64),
pointer,
dimension(:) :: wptr
2878 real(real64),
allocatable,
target,
dimension(:) :: wrk
2879 class(errors),
pointer :: errmgr
2880 type(errors),
target :: deferr
2886 if (
present(err))
then
2893 if (
size(tau) /= mn)
then
2895 call errmgr%report_error(
"lq_factor_no_pivot", &
2896 "Incorrectly sized input array TAU, argument 2.", &
2897 la_array_size_error)
2902 call dgelqf(m, n, a, m, tau, temp, -1, flag)
2903 lwork = int(temp(1), int32)
2904 if (
present(olwork))
then
2910 if (
present(work))
then
2911 if (
size(work) < lwork)
then
2913 call errmgr%report_error(
"lq_factor_no_pivot", &
2914 "Incorrectly sized input array WORK, argument 3.", &
2915 la_array_size_error)
2918 wptr => work(1:lwork)
2920 allocate(wrk(lwork), stat = istat)
2921 if (istat /= 0)
then
2923 call errmgr%report_error(
"lq_factor_no_pivot", &
2924 "Insufficient memory available.", &
2925 la_out_of_memory_error)
2932 call dgelqf(m, n, a, m, tau, wptr, lwork, flag)
2936 module subroutine lq_factor_no_pivot_cmplx(a, tau, work, olwork, err)
2938 complex(real64),
intent(inout),
dimension(:,:) :: a
2939 complex(real64),
intent(out),
dimension(:) :: tau
2940 complex(real64),
intent(out),
target,
dimension(:),
optional :: work
2941 integer(int32),
intent(out),
optional :: olwork
2942 class(errors),
intent(inout),
optional,
target :: err
2945 integer(int32) :: m, n, mn, istat, lwork, flag
2946 complex(real64),
dimension(1) :: temp
2947 complex(real64),
pointer,
dimension(:) :: wptr
2948 complex(real64),
allocatable,
target,
dimension(:) :: wrk
2949 class(errors),
pointer :: errmgr
2950 type(errors),
target :: deferr
2956 if (
present(err))
then
2963 if (
size(tau) /= mn)
then
2965 call errmgr%report_error(
"lq_factor_no_pivot_cmplx", &
2966 "Incorrectly sized input array TAU, argument 2.", &
2967 la_array_size_error)
2972 call zgelqf(m, n, a, m, tau, temp, -1, flag)
2973 lwork = int(temp(1), int32)
2974 if (
present(olwork))
then
2980 if (
present(work))
then
2981 if (
size(work) < lwork)
then
2983 call errmgr%report_error(
"lq_factor_no_pivot_cmplx", &
2984 "Incorrectly sized input array WORK, argument 3.", &
2985 la_array_size_error)
2988 wptr => work(1:lwork)
2990 allocate(wrk(lwork), stat = istat)
2991 if (istat /= 0)
then
2993 call errmgr%report_error(
"lq_factor_no_pivot_cmplx", &
2994 "Insufficient memory available.", &
2995 la_out_of_memory_error)
3002 call zgelqf(m, n, a, m, tau, wptr, lwork, flag)
3006 module subroutine form_lq_no_pivot(l, tau, q, work, olwork, err)
3008 real(real64),
intent(inout),
dimension(:,:) :: l
3009 real(real64),
intent(in),
dimension(:) :: tau
3010 real(real64),
intent(out),
dimension(:,:) :: q
3011 real(real64),
intent(out),
target,
dimension(:),
optional :: work
3012 integer(int32),
intent(out),
optional :: olwork
3013 class(errors),
intent(inout),
optional,
target :: err
3016 real(real64),
parameter :: zero = 0.0d0
3019 integer(int32) :: i, m, n, mn, k, istat, flag, lwork
3020 real(real64),
pointer,
dimension(:) :: wptr
3021 real(real64),
allocatable,
target,
dimension(:) :: wrk
3022 real(real64),
dimension(1) :: temp
3023 class(errors),
pointer :: errmgr
3024 type(errors),
target :: deferr
3025 character(len = 128) :: errmsg
3031 if (
present(err))
then
3041 else if (
size(tau) /= mn)
then
3043 else if (
size(q, 1) /= n .or.
size(q, 2) /= n)
then
3048 write(errmsg, 100)
"Input number ", flag, &
3049 " is not sized correctly."
3050 call errmgr%report_error(
"form_lq_no_pivot", trim(errmsg), &
3051 la_array_size_error)
3056 call dorglq(n, n, mn, q, n, tau, temp, -1, flag)
3057 lwork = int(temp(1), int32)
3058 if (
present(olwork))
then
3064 if (
present(work))
then
3065 if (
size(work) < lwork)
then
3067 call errmgr%report_error(
"form_lq_no_pivot", &
3068 "Incorrectly sized input array WORK, argument 4.", &
3069 la_array_size_error)
3072 wptr => work(1:lwork)
3074 allocate(wrk(lwork), stat = istat)
3075 if (istat /= 0)
then
3077 call errmgr%report_error(
"form_lq_no_pivot", &
3078 "Insufficient memory available.", &
3079 la_out_of_memory_error)
3088 q(1:j-1,j) = l(1:k,j)
3093 call dorglq(n, n, mn, q, n, tau, wptr, lwork, flag)
3100 module subroutine form_lq_no_pivot_cmplx(l, tau, q, work, olwork, err)
3102 complex(real64),
intent(inout),
dimension(:,:) :: l
3103 complex(real64),
intent(in),
dimension(:) :: tau
3104 complex(real64),
intent(out),
dimension(:,:) :: q
3105 complex(real64),
intent(out),
target,
dimension(:),
optional :: work
3106 integer(int32),
intent(out),
optional :: olwork
3107 class(errors),
intent(inout),
optional,
target :: err
3110 complex(real64),
parameter :: zero = (0.0d0, 0.0d0)
3113 integer(int32) :: i, m, n, mn, k, istat, flag, lwork
3114 complex(real64),
pointer,
dimension(:) :: wptr
3115 complex(real64),
allocatable,
target,
dimension(:) :: wrk
3116 complex(real64),
dimension(1) :: temp
3117 class(errors),
pointer :: errmgr
3118 type(errors),
target :: deferr
3119 character(len = 128) :: errmsg
3125 if (
present(err))
then
3135 else if (
size(tau) /= mn)
then
3137 else if (
size(q, 1) /= n .or.
size(q, 2) /= n)
then
3142 write(errmsg, 100)
"Input number ", flag, &
3143 " is not sized correctly."
3144 call errmgr%report_error(
"form_lq_no_pivot_cmplx", trim(errmsg), &
3145 la_array_size_error)
3150 call zunglq(n, n, mn, q, n, tau, temp, -1, flag)
3151 lwork = int(temp(1), int32)
3152 if (
present(olwork))
then
3158 if (
present(work))
then
3159 if (
size(work) < lwork)
then
3161 call errmgr%report_error(
"form_lq_no_pivot_cmplx", &
3162 "Incorrectly sized input array WORK, argument 4.", &
3163 la_array_size_error)
3166 wptr => work(1:lwork)
3168 allocate(wrk(lwork), stat = istat)
3169 if (istat /= 0)
then
3171 call errmgr%report_error(
"form_lq_no_pivot_cmplx", &
3172 "Insufficient memory available.", &
3173 la_out_of_memory_error)
3182 q(1:j-1,j) = l(1:k,j)
3187 call zunglq(n, n, mn, q, n, tau, wptr, lwork, flag)
3194 module subroutine mult_lq_mtx(lside, trans, a, tau, c, work, olwork, err)
3196 logical,
intent(in) :: lside, trans
3197 real(real64),
intent(in),
dimension(:,:) :: a
3198 real(real64),
intent(in),
dimension(:) :: tau
3199 real(real64),
intent(inout),
dimension(:,:) :: c
3200 real(real64),
intent(out),
target,
dimension(:),
optional :: work
3201 integer(int32),
intent(out),
optional :: olwork
3202 class(errors),
intent(inout),
optional,
target :: err
3205 character :: side, t
3206 integer(int32) :: m, n, k, ncola, istat, flag, lwork
3207 real(real64),
pointer,
dimension(:) :: wptr
3208 real(real64),
allocatable,
target,
dimension(:) :: wrk
3209 real(real64),
dimension(1) :: temp
3210 class(errors),
pointer :: errmgr
3211 type(errors),
target :: deferr
3212 character(len = 128) :: errmsg
3230 if (
present(err))
then
3238 if (
size(a, 1) /= k .or.
size(a, 2) /= ncola)
then
3243 write(errmsg, 100)
"Input number ", flag, &
3244 " is not sized correctly."
3245 call errmgr%report_error(
"mult_lq_mtx", trim(errmsg), &
3246 la_array_size_error)
3251 call dormlq(side, t, m, n, k, a, k, tau, c, m, temp, -1, flag)
3252 lwork = int(temp(1), int32)
3253 if (
present(olwork))
then
3259 if (
present(work))
then
3260 if (
size(work) < lwork)
then
3262 call errmgr%report_error(
"mult_lq_mtx", &
3263 "Incorrectly sized input array WORK, argument 6.", &
3264 la_array_size_error)
3267 wptr => work(1:lwork)
3269 allocate(wrk(lwork), stat = istat)
3270 if (istat /= 0)
then
3272 call errmgr%report_error(
"mult_lq_mtx", &
3273 "Insufficient memory available.", &
3274 la_out_of_memory_error)
3281 call dormlq(side, t, m, n, k, a, k, tau, c, m, wptr, lwork, flag)
3288 module subroutine mult_lq_mtx_cmplx(lside, trans, a, tau, c, work, olwork, err)
3290 logical,
intent(in) :: lside, trans
3291 complex(real64),
intent(in),
dimension(:,:) :: a
3292 complex(real64),
intent(in),
dimension(:) :: tau
3293 complex(real64),
intent(inout),
dimension(:,:) :: c
3294 complex(real64),
intent(out),
target,
dimension(:),
optional :: work
3295 integer(int32),
intent(out),
optional :: olwork
3296 class(errors),
intent(inout),
optional,
target :: err
3299 character :: side, t
3300 integer(int32) :: m, n, k, ncola, istat, flag, lwork
3301 complex(real64),
pointer,
dimension(:) :: wptr
3302 complex(real64),
allocatable,
target,
dimension(:) :: wrk
3303 complex(real64),
dimension(1) :: temp
3304 class(errors),
pointer :: errmgr
3305 type(errors),
target :: deferr
3306 character(len = 128) :: errmsg
3324 if (
present(err))
then
3332 if (
size(a, 1) /= k .or.
size(a, 2) /= ncola)
then
3337 write(errmsg, 100)
"Input number ", flag, &
3338 " is not sized correctly."
3339 call errmgr%report_error(
"mult_lq_mtx_cmplx", trim(errmsg), &
3340 la_array_size_error)
3345 call zunmlq(side, t, m, n, k, a, k, tau, c, m, temp, -1, flag)
3346 lwork = int(temp(1), int32)
3347 if (
present(olwork))
then
3353 if (
present(work))
then
3354 if (
size(work) < lwork)
then
3356 call errmgr%report_error(
"mult_lq_mtx_cmplx", &
3357 "Incorrectly sized input array WORK, argument 6.", &
3358 la_array_size_error)
3361 wptr => work(1:lwork)
3363 allocate(wrk(lwork), stat = istat)
3364 if (istat /= 0)
then
3366 call errmgr%report_error(
"mult_lq_mtx_cmplx", &
3367 "Insufficient memory available.", &
3368 la_out_of_memory_error)
3375 call zunmlq(side, t, m, n, k, a, k, tau, c, m, wptr, lwork, flag)
3382 module subroutine mult_lq_vec(trans, a, tau, c, work, olwork, err)
3384 logical,
intent(in) :: trans
3385 real(real64),
intent(in),
dimension(:,:) :: a
3386 real(real64),
intent(in),
dimension(:) :: tau
3387 real(real64),
intent(inout),
dimension(:) :: c
3388 real(real64),
intent(out),
target,
dimension(:),
optional :: work
3389 integer(int32),
intent(out),
optional :: olwork
3390 class(errors),
intent(inout),
optional,
target :: err
3393 character :: side, t
3394 integer(int32) :: m, n, k, istat, flag, lwork
3395 real(real64),
pointer,
dimension(:) :: wptr
3396 real(real64),
allocatable,
target,
dimension(:) :: wrk
3397 real(real64),
dimension(1) :: temp
3398 class(errors),
pointer :: errmgr
3399 type(errors),
target :: deferr
3400 character(len = 128) :: errmsg
3412 if (
present(err))
then
3420 if (
size(a, 1) /= k .or.
size(a, 2) /= m)
then
3425 write(errmsg, 100)
"Input number ", flag, &
3426 " is not sized correctly."
3427 call errmgr%report_error(
"mult_lq_vec", trim(errmsg), &
3428 la_array_size_error)
3433 call dormlq(side, t, m, n, k, a, k, tau, c, m, temp, -1, flag)
3434 lwork = int(temp(1), int32)
3435 if (
present(olwork))
then
3441 if (
present(work))
then
3442 if (
size(work) < lwork)
then
3444 call errmgr%report_error(
"mult_lq_vec", &
3445 "Incorrectly sized input array WORK, argument 6.", &
3446 la_array_size_error)
3449 wptr => work(1:lwork)
3451 allocate(wrk(lwork), stat = istat)
3452 if (istat /= 0)
then
3454 call errmgr%report_error(
"mult_lq_vec", &
3455 "Insufficient memory available.", &
3456 la_out_of_memory_error)
3463 call dormlq(side, t, m, n, k, a, k, tau, c, m, wptr, lwork, flag)
3470 module subroutine mult_lq_vec_cmplx(trans, a, tau, c, work, olwork, err)
3472 logical,
intent(in) :: trans
3473 complex(real64),
intent(in),
dimension(:,:) :: a
3474 complex(real64),
intent(in),
dimension(:) :: tau
3475 complex(real64),
intent(inout),
dimension(:) :: c
3476 complex(real64),
intent(out),
target,
dimension(:),
optional :: work
3477 integer(int32),
intent(out),
optional :: olwork
3478 class(errors),
intent(inout),
optional,
target :: err
3481 character :: side, t
3482 integer(int32) :: m, n, k, istat, flag, lwork
3483 complex(real64),
pointer,
dimension(:) :: wptr
3484 complex(real64),
allocatable,
target,
dimension(:) :: wrk
3485 complex(real64),
dimension(1) :: temp
3486 class(errors),
pointer :: errmgr
3487 type(errors),
target :: deferr
3488 character(len = 128) :: errmsg
3500 if (
present(err))
then
3508 if (
size(a, 1) /= k .or.
size(a, 2) /= m)
then
3513 write(errmsg, 100)
"Input number ", flag, &
3514 " is not sized correctly."
3515 call errmgr%report_error(
"mult_lq_vec_cmplx", trim(errmsg), &
3516 la_array_size_error)
3521 call zunmlq(side, t, m, n, k, a, k, tau, c, m, temp, -1, flag)
3522 lwork = int(temp(1), int32)
3523 if (
present(olwork))
then
3529 if (
present(work))
then
3530 if (
size(work) < lwork)
then
3532 call errmgr%report_error(
"mult_lq_vec_cmplx", &
3533 "Incorrectly sized input array WORK, argument 6.", &
3534 la_array_size_error)
3537 wptr => work(1:lwork)
3539 allocate(wrk(lwork), stat = istat)
3540 if (istat /= 0)
then
3542 call errmgr%report_error(
"mult_lq_vec_cmplx", &
3543 "Insufficient memory available.", &
3544 la_out_of_memory_error)
3551 call zunmlq(side, t, m, n, k, a, k, tau, c, m, wptr, lwork, flag)
Provides a set of common linear algebra routines.