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 = :),
allocatable :: 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)
52 allocate(
character(len = 256) :: errmsg)
54 "Singular matrix encountered (row ", flag,
")"
55 call errmgr%report_warning(
"lu_factor_dbl", trim(errmsg), &
56 la_singular_matrix_error)
64 module subroutine lu_factor_cmplx(a, ipvt, err)
66 complex(real64),
intent(inout),
dimension(:,:) :: a
67 integer(int32),
intent(out),
dimension(:) :: ipvt
68 class(errors),
intent(inout),
optional,
target :: err
71 integer(int32) :: m, n, mn, flag
72 class(errors),
pointer :: errmgr
73 type(errors),
target :: deferr
74 character(len = :),
allocatable :: errmsg
80 if (
present(err))
then
88 if (
size(ipvt) /= mn)
then
90 call errmgr%report_error(
"lu_factor_cmplx", &
91 "Incorrectly sized input array IPVT, argument 2.", &
97 call zgetrf(m, n, a, m, ipvt, flag)
104 allocate(
character(len = 256) :: errmsg)
106 "Singular matrix encountered (row ", flag,
")"
107 call errmgr%report_warning(
"lu_factor_cmplx", trim(errmsg), &
108 la_singular_matrix_error)
116 module subroutine form_lu_all(lu, ipvt, u, p, err)
118 real(real64),
intent(inout),
dimension(:,:) :: lu
119 integer(int32),
intent(in),
dimension(:) :: ipvt
120 real(real64),
intent(out),
dimension(:,:) :: u, p
121 class(errors),
intent(inout),
optional,
target :: err
124 integer(int32) :: j, jp, n, flag
125 class(errors),
pointer :: errmgr
126 type(errors),
target :: deferr
127 character(len = :),
allocatable :: errmsg
130 real(real64),
parameter :: zero = 0.0d0
131 real(real64),
parameter :: one = 1.0d0
135 if (
present(err))
then
143 if (
size(lu, 2) /= n)
then
145 else if (
size(ipvt) /= n)
then
147 else if (
size(u, 1) /= n .or.
size(u, 2) /= n)
then
149 else if (
size(p, 1) /= n .or.
size(p, 2) /= n)
then
154 allocate(
character(len = 256) :: errmsg)
155 write(errmsg, 100)
"Input number ", flag, &
156 " is not sized correctly."
157 call errmgr%report_error(
"form_lu_all", trim(errmsg), &
163 call dlaset(
'A', n, n, zero, one, p, n)
169 if (j /= jp)
call swap(p(j,1:n), p(jp,1:n))
175 if (j > 1) lu(1:j-1,j) = zero
184 module subroutine form_lu_all_cmplx(lu, ipvt, u, p, err)
186 complex(real64),
intent(inout),
dimension(:,:) :: lu
187 integer(int32),
intent(in),
dimension(:) :: ipvt
188 complex(real64),
intent(out),
dimension(:,:) :: u
189 real(real64),
intent(out),
dimension(:,:) :: p
190 class(errors),
intent(inout),
optional,
target :: err
193 integer(int32) :: j, jp, n, flag
194 class(errors),
pointer :: errmgr
195 type(errors),
target :: deferr
196 character(len = :),
allocatable :: errmsg
199 real(real64),
parameter :: zero = 0.0d0
200 real(real64),
parameter :: one = 1.0d0
201 complex(real64),
parameter :: c_zero = (0.0d0, 0.0d0)
202 complex(real64),
parameter :: c_one = (1.0d0, 0.0d0)
206 if (
present(err))
then
214 if (
size(lu, 2) /= n)
then
216 else if (
size(ipvt) /= n)
then
218 else if (
size(u, 1) /= n .or.
size(u, 2) /= n)
then
220 else if (
size(p, 1) /= n .or.
size(p, 2) /= n)
then
225 allocate(
character(len = 256) :: errmsg)
226 write(errmsg, 100)
"Input number ", flag, &
227 " is not sized correctly."
228 call errmgr%report_error(
"form_lu_all_cmplx", trim(errmsg), &
234 call dlaset(
'A', n, n, zero, one, p, n)
240 if (j /= jp)
call swap(p(j,1:n), p(jp,1:n))
246 if (j > 1) lu(1:j-1,j) = c_zero
255 module subroutine form_lu_only(lu, u, err)
257 real(real64),
intent(inout),
dimension(:,:) :: lu
258 real(real64),
intent(out),
dimension(:,:) :: u
259 class(errors),
intent(inout),
optional,
target :: err
262 integer(int32) :: j, n, flag
263 class(errors),
pointer :: errmgr
264 type(errors),
target :: deferr
265 character(len = :),
allocatable :: errmsg
268 real(real64),
parameter :: zero = 0.0d0
269 real(real64),
parameter :: one = 1.0d0
273 if (
present(err))
then
281 if (
size(lu, 2) /= n)
then
283 else if (
size(u, 1) /= n .or.
size(u, 2) /= n)
then
288 allocate(
character(len = 256) :: errmsg)
289 write(errmsg, 100)
"Input number ", flag, &
290 " is not sized correctly."
291 call errmgr%report_error(
"form_lu_only", trim(errmsg), &
302 if (j > 1) lu(1:j-1,j) = zero
311 module subroutine form_lu_only_cmplx(lu, u, err)
313 complex(real64),
intent(inout),
dimension(:,:) :: lu
314 complex(real64),
intent(out),
dimension(:,:) :: u
315 class(errors),
intent(inout),
optional,
target :: err
318 integer(int32) :: j, n, flag
319 class(errors),
pointer :: errmgr
320 type(errors),
target :: deferr
321 character(len = :),
allocatable :: errmsg
324 complex(real64),
parameter :: zero = (0.0d0, 0.0d0)
325 complex(real64),
parameter :: one = (1.0d0, 0.0d0)
329 if (
present(err))
then
337 if (
size(lu, 2) /= n)
then
339 else if (
size(u, 1) /= n .or.
size(u, 2) /= n)
then
344 allocate(
character(len = 256) :: errmsg)
345 write(errmsg, 100)
"Input number ", flag, &
346 " is not sized correctly."
347 call errmgr%report_error(
"form_lu_only_cmplx", trim(errmsg), &
358 if (j > 1) lu(1:j-1,j) = zero
369 module subroutine qr_factor_no_pivot(a, tau, work, olwork, err)
371 real(real64),
intent(inout),
dimension(:,:) :: a
372 real(real64),
intent(out),
dimension(:) :: tau
373 real(real64),
intent(out),
target,
dimension(:),
optional :: work
374 integer(int32),
intent(out),
optional :: olwork
375 class(errors),
intent(inout),
optional,
target :: err
378 integer(int32) :: m, n, mn, istat, lwork, flag
379 real(real64),
dimension(1) :: temp
380 real(real64),
pointer,
dimension(:) :: wptr
381 real(real64),
allocatable,
target,
dimension(:) :: wrk
382 class(errors),
pointer :: errmgr
383 type(errors),
target :: deferr
389 if (
present(err))
then
396 if (
size(tau) /= mn)
then
398 call errmgr%report_error(
"qr_factor_no_pivot", &
399 "Incorrectly sized input array TAU, argument 2.", &
405 call dgeqrf(m, n, a, m, tau, temp, -1, flag)
406 lwork = int(temp(1), int32)
407 if (
present(olwork))
then
413 if (
present(work))
then
414 if (
size(work) < lwork)
then
416 call errmgr%report_error(
"qr_factor_no_pivot", &
417 "Incorrectly sized input array WORK, argument 3.", &
421 wptr => work(1:lwork)
423 allocate(wrk(lwork), stat = istat)
426 call errmgr%report_error(
"qr_factor_no_pivot", &
427 "Insufficient memory available.", &
428 la_out_of_memory_error)
435 call dgeqrf(m, n, a, m, tau, wptr, lwork, flag)
439 module subroutine qr_factor_no_pivot_cmplx(a, tau, work, olwork, err)
441 complex(real64),
intent(inout),
dimension(:,:) :: a
442 complex(real64),
intent(out),
dimension(:) :: tau
443 complex(real64),
intent(out),
target,
dimension(:),
optional :: work
444 integer(int32),
intent(out),
optional :: olwork
445 class(errors),
intent(inout),
optional,
target :: err
448 integer(int32) :: m, n, mn, istat, lwork, flag
449 complex(real64),
dimension(1) :: temp
450 complex(real64),
pointer,
dimension(:) :: wptr
451 complex(real64),
allocatable,
target,
dimension(:) :: wrk
452 class(errors),
pointer :: errmgr
453 type(errors),
target :: deferr
459 if (
present(err))
then
466 if (
size(tau) /= mn)
then
468 call errmgr%report_error(
"qr_factor_no_pivot_cmplx", &
469 "Incorrectly sized input array TAU, argument 2.", &
475 call zgeqrf(m, n, a, m, tau, temp, -1, flag)
476 lwork = int(temp(1), int32)
477 if (
present(olwork))
then
483 if (
present(work))
then
484 if (
size(work) < lwork)
then
486 call errmgr%report_error(
"qr_factor_no_pivot_cmplx", &
487 "Incorrectly sized input array WORK, argument 3.", &
491 wptr => work(1:lwork)
493 allocate(wrk(lwork), stat = istat)
496 call errmgr%report_error(
"qr_factor_no_pivot_cmplx", &
497 "Insufficient memory available.", &
498 la_out_of_memory_error)
505 call zgeqrf(m, n, a, m, tau, wptr, lwork, flag)
509 module subroutine qr_factor_pivot(a, tau, jpvt, work, olwork, err)
511 real(real64),
intent(inout),
dimension(:,:) :: a
512 real(real64),
intent(out),
dimension(:) :: tau
513 integer(int32),
intent(inout),
dimension(:) :: jpvt
514 real(real64),
intent(out),
target,
dimension(:),
optional :: work
515 integer(int32),
intent(out),
optional :: olwork
516 class(errors),
intent(inout),
optional,
target :: err
519 integer(int32) :: m, n, mn, istat, lwork, flag
520 real(real64),
dimension(1) :: temp
521 real(real64),
pointer,
dimension(:) :: wptr
522 real(real64),
allocatable,
target,
dimension(:) :: wrk
523 class(errors),
pointer :: errmgr
524 type(errors),
target :: deferr
525 character(len = :),
allocatable :: errmsg
531 if (
present(err))
then
539 if (
size(tau) /= mn)
then
541 else if (
size(jpvt) /= n)
then
546 allocate(
character(len = 256) :: errmsg)
547 write(errmsg, 100)
"Input number ", flag, &
548 " is not sized correctly."
549 call errmgr%report_error(
"qr_factor_pivot", trim(errmsg), &
555 call dgeqp3(m, n, a, m, jpvt, tau, temp, -1, flag)
556 lwork = int(temp(1), int32)
557 if (
present(olwork))
then
563 if (
present(work))
then
564 if (
size(work) < lwork)
then
566 call errmgr%report_error(
"qr_factor_pivot", &
567 "Incorrectly sized input array WORK, argument 4.", &
571 wptr => work(1:lwork)
573 allocate(wrk(lwork), stat = istat)
576 call errmgr%report_error(
"qr_factor_pivot", &
577 "Insufficient memory available.", &
578 la_out_of_memory_error)
585 call dgeqp3(m, n, a, m, jpvt, tau, wptr, lwork, flag)
588 if (
allocated(wrk))
deallocate(wrk)
595 module subroutine qr_factor_pivot_cmplx(a, tau, jpvt, work, olwork, rwork, &
598 complex(real64),
intent(inout),
dimension(:,:) :: a
599 complex(real64),
intent(out),
dimension(:) :: tau
600 integer(int32),
intent(inout),
dimension(:) :: jpvt
601 complex(real64),
intent(out),
target,
dimension(:),
optional :: work
602 integer(int32),
intent(out),
optional :: olwork
603 real(real64),
intent(out),
target,
dimension(:),
optional :: rwork
604 class(errors),
intent(inout),
optional,
target :: err
607 integer(int32) :: m, n, mn, istat, lwork, flag
608 complex(real64),
dimension(1) :: temp
609 complex(real64),
pointer,
dimension(:) :: wptr
610 complex(real64),
allocatable,
target,
dimension(:) :: wrk
611 real(real64),
pointer,
dimension(:) :: rptr
612 real(real64),
allocatable,
target,
dimension(:) :: rwrk
613 class(errors),
pointer :: errmgr
614 type(errors),
target :: deferr
615 character(len = :),
allocatable :: errmsg
621 if (
present(err))
then
629 if (
size(tau) /= mn)
then
631 else if (
size(jpvt) /= n)
then
636 allocate(
character(len = 256) :: errmsg)
637 write(errmsg, 100)
"Input number ", flag, &
638 " is not sized correctly."
639 call errmgr%report_error(
"qr_factor_pivot_cmplx", trim(errmsg), &
643 if (
present(rwork))
then
644 if (
size(rwork) < 2 * n)
then
645 call errmgr%report_error(
"qr_factor_pivot_cmplx", &
646 "Incorrectly sized input array RWORK, argument 6.", &
652 allocate(rwrk(2 * n), stat = flag)
654 call errmgr%report_error(
"qr_factor_pivot_cmplx", &
655 "Insufficient memory available.", &
656 la_out_of_memory_error)
663 call zgeqp3(m, n, a, m, jpvt, tau, temp, -1, rptr, flag)
664 lwork = int(temp(1), int32)
665 if (
present(olwork))
then
671 if (
present(work))
then
672 if (
size(work) < lwork)
then
674 call errmgr%report_error(
"qr_factor_pivot_cmplx", &
675 "Incorrectly sized input array WORK, argument 4.", &
679 wptr => work(1:lwork)
681 allocate(wrk(lwork), stat = istat)
684 call errmgr%report_error(
"qr_factor_pivot_cmplx", &
685 "Insufficient memory available.", &
686 la_out_of_memory_error)
693 call zgeqp3(m, n, a, m, jpvt, tau, wptr, lwork, rptr, flag)
696 if (
allocated(wrk))
deallocate(wrk)
703 module subroutine form_qr_no_pivot(r, tau, q, work, olwork, err)
705 real(real64),
intent(inout),
dimension(:,:) :: r
706 real(real64),
intent(in),
dimension(:) :: tau
707 real(real64),
intent(out),
dimension(:,:) :: q
708 real(real64),
intent(out),
target,
dimension(:),
optional :: work
709 integer(int32),
intent(out),
optional :: olwork
710 class(errors),
intent(inout),
optional,
target :: err
713 real(real64),
parameter :: zero = 0.0d0
716 integer(int32) :: j, m, n, mn, qcol, istat, flag, lwork
717 real(real64),
pointer,
dimension(:) :: wptr
718 real(real64),
allocatable,
target,
dimension(:) :: wrk
719 real(real64),
dimension(1) :: temp
720 class(errors),
pointer :: errmgr
721 type(errors),
target :: deferr
722 character(len = :),
allocatable :: errmsg
729 if (
present(err))
then
737 if (
size(tau) /= mn)
then
739 else if (
size(q, 1) /= m .or. (qcol /= m .and. qcol /= n))
then
741 else if (qcol == n .and. m < n)
then
746 allocate(
character(len = 256) :: errmsg)
747 write(errmsg, 100)
"Input number ", flag, &
748 " is not sized correctly."
749 call errmgr%report_error(
"form_qr_no_pivot", trim(errmsg), &
755 call dorgqr(m, qcol, mn, q, m, tau, temp, -1, flag)
756 lwork = int(temp(1), int32)
757 if (
present(olwork))
then
763 if (
present(work))
then
764 if (
size(work) < lwork)
then
766 call errmgr%report_error(
"form_qr_no_pivot", &
767 "Incorrectly sized input array WORK, argument 4.", &
771 wptr => work(1:lwork)
773 allocate(wrk(lwork), stat = istat)
776 call errmgr%report_error(
"form_qr_no_pivot", &
777 "Insufficient memory available.", &
778 la_out_of_memory_error)
787 q(j+1:m,j) = r(j+1:m,j)
792 call dorgqr(m, qcol, mn, q, m, tau, wptr, lwork, flag)
795 if (
allocated(wrk))
deallocate(wrk)
802 module subroutine form_qr_no_pivot_cmplx(r, tau, q, work, olwork, err)
804 complex(real64),
intent(inout),
dimension(:,:) :: r
805 complex(real64),
intent(in),
dimension(:) :: tau
806 complex(real64),
intent(out),
dimension(:,:) :: q
807 complex(real64),
intent(out),
target,
dimension(:),
optional :: work
808 integer(int32),
intent(out),
optional :: olwork
809 class(errors),
intent(inout),
optional,
target :: err
812 complex(real64),
parameter :: zero = (0.0d0, 0.0d0)
815 integer(int32) :: j, m, n, mn, qcol, istat, flag, lwork
816 complex(real64),
pointer,
dimension(:) :: wptr
817 complex(real64),
allocatable,
target,
dimension(:) :: wrk
818 complex(real64),
dimension(1) :: temp
819 class(errors),
pointer :: errmgr
820 type(errors),
target :: deferr
821 character(len = :),
allocatable :: errmsg
828 if (
present(err))
then
836 if (
size(tau) /= mn)
then
838 else if (
size(q, 1) /= m .or. (qcol /= m .and. qcol /= n))
then
840 else if (qcol == n .and. m < n)
then
845 allocate(
character(len = 256) :: errmsg)
846 write(errmsg, 100)
"Input number ", flag, &
847 " is not sized correctly."
848 call errmgr%report_error(
"form_qr_no_pivot_cmplx", trim(errmsg), &
854 call zungqr(m, qcol, mn, q, m, tau, temp, -1, flag)
855 lwork = int(temp(1), int32)
856 if (
present(olwork))
then
862 if (
present(work))
then
863 if (
size(work) < lwork)
then
865 call errmgr%report_error(
"form_qr_no_pivot_cmplx", &
866 "Incorrectly sized input array WORK, argument 4.", &
870 wptr => work(1:lwork)
872 allocate(wrk(lwork), stat = istat)
875 call errmgr%report_error(
"form_qr_no_pivot_cmplx", &
876 "Insufficient memory available.", &
877 la_out_of_memory_error)
886 q(j+1:m,j) = r(j+1:m,j)
891 call zungqr(m, qcol, mn, q, m, tau, wptr, lwork, flag)
894 if (
allocated(wrk))
deallocate(wrk)
901 module subroutine form_qr_pivot(r, tau, pvt, q, p, work, olwork, err)
903 real(real64),
intent(inout),
dimension(:,:) :: r
904 real(real64),
intent(in),
dimension(:) :: tau
905 integer(int32),
intent(in),
dimension(:) :: pvt
906 real(real64),
intent(out),
dimension(:,:) :: q, p
907 real(real64),
intent(out),
target,
dimension(:),
optional :: work
908 integer(int32),
intent(out),
optional :: olwork
909 class(errors),
intent(inout),
optional,
target :: err
912 real(real64),
parameter :: zero = 0.0d0
913 real(real64),
parameter :: one = 1.0d0
916 integer(int32) :: j, jp, m, n, mn, flag
917 class(errors),
pointer :: errmgr
918 type(errors),
target :: deferr
919 character(len = :),
allocatable :: errmsg
925 if (
present(err))
then
933 if (
size(tau) /= mn)
then
935 else if (
size(pvt) /= n)
then
937 else if (
size(q, 1) /= m .or. &
938 (
size(q, 2) /= m .and.
size(q, 2) /= n))
then
940 else if (
size(q, 2) == n .and. m < n)
then
942 else if (
size(p, 1) /= n .or.
size(p, 2) /= n)
then
947 allocate(
character(len = 256) :: errmsg)
948 write(errmsg, 100)
"Input number ", flag, &
949 " is not sized correctly."
950 call errmgr%report_error(
"form_qr_pivot", trim(errmsg), &
956 call form_qr_no_pivot(r, tau, q, work = work, olwork = olwork, &
958 if (
present(olwork))
return
959 if (errmgr%has_error_occurred())
return
973 module subroutine form_qr_pivot_cmplx(r, tau, pvt, q, p, work, olwork, err)
975 complex(real64),
intent(inout),
dimension(:,:) :: r
976 complex(real64),
intent(in),
dimension(:) :: tau
977 integer(int32),
intent(in),
dimension(:) :: pvt
978 complex(real64),
intent(out),
dimension(:,:) :: q, p
979 complex(real64),
intent(out),
target,
dimension(:),
optional :: work
980 integer(int32),
intent(out),
optional :: olwork
981 class(errors),
intent(inout),
optional,
target :: err
984 complex(real64),
parameter :: zero = (0.0d0, 0.0d0)
985 complex(real64),
parameter :: one = (1.0d0, 0.0d0)
988 integer(int32) :: j, jp, m, n, mn, flag
989 class(errors),
pointer :: errmgr
990 type(errors),
target :: deferr
991 character(len = :),
allocatable :: errmsg
997 if (
present(err))
then
1005 if (
size(tau) /= mn)
then
1007 else if (
size(pvt) /= n)
then
1009 else if (
size(q, 1) /= m .or. &
1010 (
size(q, 2) /= m .and.
size(q, 2) /= n))
then
1012 else if (
size(q, 2) == n .and. m < n)
then
1014 else if (
size(p, 1) /= n .or.
size(p, 2) /= n)
then
1019 allocate(
character(len = 256) :: errmsg)
1020 write(errmsg, 100)
"Input number ", flag, &
1021 " is not sized correctly."
1022 call errmgr%report_error(
"form_qr_pivot_cmplx", trim(errmsg), &
1023 la_array_size_error)
1028 call form_qr_no_pivot_cmplx(r, tau, q, work = work, olwork = olwork, &
1030 if (
present(olwork))
return
1031 if (errmgr%has_error_occurred())
return
1045 module subroutine mult_qr_mtx(lside, trans, a, tau, c, work, olwork, err)
1047 logical,
intent(in) :: lside, trans
1048 real(real64),
intent(in),
dimension(:) :: tau
1049 real(real64),
intent(inout),
dimension(:,:) :: a, c
1050 real(real64),
intent(out),
target,
dimension(:),
optional :: work
1051 integer(int32),
intent(out),
optional :: olwork
1052 class(errors),
intent(inout),
optional,
target :: err
1055 real(real64),
parameter :: one = 1.0d0
1058 character :: side, t
1059 integer(int32) :: m, n, k, nrowa, istat, flag, lwork
1060 real(real64),
pointer,
dimension(:) :: wptr
1061 real(real64),
allocatable,
target,
dimension(:) :: wrk
1062 real(real64),
dimension(1) :: temp
1063 class(errors),
pointer :: errmgr
1064 type(errors),
target :: deferr
1065 character(len = :),
allocatable :: errmsg
1083 if (
present(err))
then
1093 if (
size(a, 1) /= m .or.
size(a, 2) < k) flag = 3
1096 if (
size(a, 1) /= n .or.
size(a, 2) < k) flag = 3
1100 allocate(
character(len = 256) :: errmsg)
1101 write(errmsg, 100)
"Input number ", flag, &
1102 " is not sized correctly."
1103 call errmgr%report_error(
"mult_qr_mtx", trim(errmsg), &
1104 la_array_size_error)
1109 call dormqr(side, t, m, n, k, a, nrowa, tau, c, m, temp, -1, flag)
1110 lwork = int(temp(1), int32)
1111 if (
present(olwork))
then
1117 if (
present(work))
then
1118 if (
size(work) < lwork)
then
1120 call errmgr%report_error(
"mult_qr_mtx", &
1121 "Incorrectly sized input array WORK, argument 6.", &
1122 la_array_size_error)
1125 wptr => work(1:lwork)
1127 allocate(wrk(lwork), stat = istat)
1128 if (istat /= 0)
then
1130 call errmgr%report_error(
"mult_qr_mtx", &
1131 "Insufficient memory available.", &
1132 la_out_of_memory_error)
1139 call dormqr(side, t, m, n, k, a, nrowa, tau, c, m, wptr, lwork, flag)
1146 module subroutine mult_qr_mtx_cmplx(lside, trans, a, tau, c, work, olwork, err)
1148 logical,
intent(in) :: lside, trans
1149 complex(real64),
intent(in),
dimension(:) :: tau
1150 complex(real64),
intent(inout),
dimension(:,:) :: a, c
1151 complex(real64),
intent(out),
target,
dimension(:),
optional :: work
1152 integer(int32),
intent(out),
optional :: olwork
1153 class(errors),
intent(inout),
optional,
target :: err
1156 complex(real64),
parameter :: one = (1.0d0, 0.0d0)
1159 character :: side, t
1160 integer(int32) :: m, n, k, nrowa, istat, flag, lwork
1161 complex(real64),
pointer,
dimension(:) :: wptr
1162 complex(real64),
allocatable,
target,
dimension(:) :: wrk
1163 complex(real64),
dimension(1) :: temp
1164 class(errors),
pointer :: errmgr
1165 type(errors),
target :: deferr
1166 character(len = :),
allocatable :: errmsg
1184 if (
present(err))
then
1194 if (
size(a, 1) /= m .or.
size(a, 2) < k) flag = 3
1197 if (
size(a, 1) /= n .or.
size(a, 2) < k) flag = 3
1201 allocate(
character(len = 256) :: errmsg)
1202 write(errmsg, 100)
"Input number ", flag, &
1203 " is not sized correctly."
1204 call errmgr%report_error(
"mult_qr_mtx_cmplx", trim(errmsg), &
1205 la_array_size_error)
1210 call zunmqr(side, t, m, n, k, a, nrowa, tau, c, m, temp, -1, flag)
1211 lwork = int(temp(1), int32)
1212 if (
present(olwork))
then
1218 if (
present(work))
then
1219 if (
size(work) < lwork)
then
1221 call errmgr%report_error(
"mult_qr_mtx_cmplx", &
1222 "Incorrectly sized input array WORK, argument 6.", &
1223 la_array_size_error)
1226 wptr => work(1:lwork)
1228 allocate(wrk(lwork), stat = istat)
1229 if (istat /= 0)
then
1231 call errmgr%report_error(
"mult_qr_mtx_cmplx", &
1232 "Insufficient memory available.", &
1233 la_out_of_memory_error)
1240 call zunmqr(side, t, m, n, k, a, nrowa, tau, c, m, wptr, lwork, flag)
1247 module subroutine mult_qr_vec(trans, a, tau, c, work, olwork, err)
1249 logical,
intent(in) :: trans
1250 real(real64),
intent(inout),
dimension(:,:) :: a
1251 real(real64),
intent(in),
dimension(:) :: tau
1252 real(real64),
intent(inout),
dimension(:) :: c
1253 real(real64),
intent(out),
target,
dimension(:),
optional :: work
1254 integer(int32),
intent(out),
optional :: olwork
1255 class(errors),
intent(inout),
optional,
target :: err
1258 real(real64),
parameter :: one = 1.0d0
1261 character :: side, t
1262 integer(int32) :: m, k, nrowa, istat, flag, lwork
1263 real(real64),
pointer,
dimension(:) :: wptr
1264 real(real64),
allocatable,
target,
dimension(:) :: wrk
1265 real(real64),
dimension(1) :: temp
1266 class(errors),
pointer :: errmgr
1267 type(errors),
target :: deferr
1268 character(len = :),
allocatable :: errmsg
1280 if (
present(err))
then
1288 if (
size(a, 1) /= m .or.
size(a, 2) < k) flag = 3
1291 allocate(
character(len = 256) :: errmsg)
1292 write(errmsg, 100)
"Input number ", flag, &
1293 " is not sized correctly."
1294 call errmgr%report_error(
"mult_qr_vec", trim(errmsg), &
1295 la_array_size_error)
1300 call dormqr(side, t, m, 1, k, a, nrowa, tau, c, m, temp, -1, flag)
1301 lwork = int(temp(1), int32)
1302 if (
present(olwork))
then
1308 if (
present(work))
then
1309 if (
size(work) < lwork)
then
1311 call errmgr%report_error(
"mult_qr_vec", &
1312 "Incorrectly sized input array WORK, argument 6.", &
1313 la_array_size_error)
1316 wptr => work(1:lwork)
1318 allocate(wrk(lwork), stat = istat)
1319 if (istat /= 0)
then
1321 call errmgr%report_error(
"mult_qr_vec", &
1322 "Insufficient memory available.", &
1323 la_out_of_memory_error)
1330 call dormqr(side, t, m, 1, k, a, nrowa, tau, c, m, wptr, lwork, flag)
1337 module subroutine mult_qr_vec_cmplx(trans, a, tau, c, work, olwork, err)
1339 logical,
intent(in) :: trans
1340 complex(real64),
intent(inout),
dimension(:,:) :: a
1341 complex(real64),
intent(in),
dimension(:) :: tau
1342 complex(real64),
intent(inout),
dimension(:) :: c
1343 complex(real64),
intent(out),
target,
dimension(:),
optional :: work
1344 integer(int32),
intent(out),
optional :: olwork
1345 class(errors),
intent(inout),
optional,
target :: err
1348 complex(real64),
parameter :: one = (1.0d0, 0.0d0)
1351 character :: side, t
1352 integer(int32) :: m, k, nrowa, istat, flag, lwork
1353 complex(real64),
pointer,
dimension(:) :: wptr
1354 complex(real64),
allocatable,
target,
dimension(:) :: wrk
1355 complex(real64),
dimension(1) :: temp
1356 class(errors),
pointer :: errmgr
1357 type(errors),
target :: deferr
1358 character(len = :),
allocatable :: errmsg
1370 if (
present(err))
then
1378 if (
size(a, 1) /= m .or.
size(a, 2) < k) flag = 3
1381 allocate(
character(len = 256) :: errmsg)
1382 write(errmsg, 100)
"Input number ", flag, &
1383 " is not sized correctly."
1384 call errmgr%report_error(
"mult_qr_vec", trim(errmsg), &
1385 la_array_size_error)
1390 call zunmqr(side, t, m, 1, k, a, nrowa, tau, c, m, temp, -1, flag)
1391 lwork = int(temp(1), int32)
1392 if (
present(olwork))
then
1398 if (
present(work))
then
1399 if (
size(work) < lwork)
then
1401 call errmgr%report_error(
"mult_qr_vec", &
1402 "Incorrectly sized input array WORK, argument 6.", &
1403 la_array_size_error)
1406 wptr => work(1:lwork)
1408 allocate(wrk(lwork), stat = istat)
1409 if (istat /= 0)
then
1411 call errmgr%report_error(
"mult_qr_vec", &
1412 "Insufficient memory available.", &
1413 la_out_of_memory_error)
1420 call zunmqr(side, t, m, 1, k, a, nrowa, tau, c, m, wptr, lwork, flag)
1427 module subroutine qr_rank1_update_dbl(q, r, u, v, work, err)
1429 real(real64),
intent(inout),
dimension(:,:) :: q, r
1430 real(real64),
intent(inout),
dimension(:) :: u, v
1431 real(real64),
intent(out),
target,
optional,
dimension(:) :: work
1432 class(errors),
intent(inout),
optional,
target :: err
1436 integer(int32) :: m, n, k, lwork, istat, flag
1437 real(real64),
pointer,
dimension(:) :: wptr
1438 real(real64),
allocatable,
target,
dimension(:) :: wrk
1439 class(errors),
pointer :: errmgr
1440 type(errors),
target :: deferr
1441 character(len = :),
allocatable :: errmsg
1447 full =
size(q, 2) == m
1449 if (
present(err))
then
1459 else if (.not.full .and.
size(q, 2) /= k)
then
1461 else if (
size(r, 1) /= m)
then
1463 else if (
size(u) /= m)
then
1465 else if (
size(v) /= n)
then
1470 allocate(
character(len = 256) :: errmsg)
1471 write(errmsg, 100)
"Input number ", flag, &
1472 " is not sized correctly."
1473 call errmgr%report_error(
"qr_rank1_update", trim(errmsg), &
1474 la_array_size_error)
1479 if (
present(work))
then
1480 if (
size(work) < lwork)
then
1482 call errmgr%report_error(
"qr_rank1_update", &
1483 "Incorrectly sized input array WORK, argument 5.", &
1484 la_array_size_error)
1487 wptr => work(1:lwork)
1489 allocate(wrk(lwork), stat = istat)
1490 if (istat /= 0)
then
1492 call errmgr%report_error(
"qr_rank1_update", &
1493 "Insufficient memory available.", &
1494 la_out_of_memory_error)
1501 call dqr1up(m, n, k, q, m, r, m, u, v, wptr)
1508 module subroutine qr_rank1_update_cmplx(q, r, u, v, work, rwork, err)
1510 complex(real64),
intent(inout),
dimension(:,:) :: q, r
1511 complex(real64),
intent(inout),
dimension(:) :: u, v
1512 complex(real64),
intent(out),
target,
optional,
dimension(:) :: work
1513 real(real64),
intent(out),
target,
optional,
dimension(:) :: rwork
1514 class(errors),
intent(inout),
optional,
target :: err
1518 integer(int32) :: m, n, k, lwork, istat, flag, lrwork
1519 complex(real64),
pointer,
dimension(:) :: wptr
1520 complex(real64),
allocatable,
target,
dimension(:) :: wrk
1521 real(real64),
pointer,
dimension(:) :: rwptr
1522 real(real64),
allocatable,
target,
dimension(:) :: rwrk
1523 class(errors),
pointer :: errmgr
1524 type(errors),
target :: deferr
1525 character(len = :),
allocatable :: errmsg
1531 full =
size(q, 2) == m
1534 if (
present(err))
then
1544 else if (.not.full .and.
size(q, 2) /= k)
then
1546 else if (
size(r, 1) /= m)
then
1548 else if (
size(u) /= m)
then
1550 else if (
size(v) /= n)
then
1555 allocate(
character(len = 256) :: errmsg)
1556 write(errmsg, 100)
"Input number ", flag, &
1557 " is not sized correctly."
1558 call errmgr%report_error(
"qr_rank1_update_cmplx", trim(errmsg), &
1559 la_array_size_error)
1564 if (
present(work))
then
1565 if (
size(work) < lwork)
then
1567 call errmgr%report_error(
"qr_rank1_update_cmplx", &
1568 "Incorrectly sized input array WORK, argument 5.", &
1569 la_array_size_error)
1572 wptr => work(1:lwork)
1574 allocate(wrk(lwork), stat = istat)
1575 if (istat /= 0)
then
1577 call errmgr%report_error(
"qr_rank1_update_cmplx", &
1578 "Insufficient memory available.", &
1579 la_out_of_memory_error)
1585 if (
present(rwork))
then
1586 if (
size(rwork) < lrwork)
then
1588 call errmgr%report_error(
"qr_rank1_update_cmplx", &
1589 "Incorrectly sized input array RWORK, argument 6.", &
1590 la_array_size_error)
1593 wptr => work(1:lrwork)
1595 allocate(rwrk(lrwork), stat = istat)
1596 if (istat /= 0)
then
1598 call errmgr%report_error(
"qr_rank1_update_cmplx", &
1599 "Insufficient memory available.", &
1600 la_out_of_memory_error)
1607 call zqr1up(m, n, k, q, m, r, m, u, v, wptr, rwptr)
1616 module subroutine cholesky_factor_dbl(a, upper, err)
1618 real(real64),
intent(inout),
dimension(:,:) :: a
1619 logical,
intent(in),
optional :: upper
1620 class(errors),
intent(inout),
optional,
target :: err
1623 real(real64),
parameter :: zero = 0.0d0
1627 integer(int32) :: i, n, flag
1628 class(errors),
pointer :: errmgr
1629 type(errors),
target :: deferr
1630 character(len = :),
allocatable :: errmsg
1634 if (
present(upper))
then
1643 if (
present(err))
then
1650 if (
size(a, 2) /= n)
then
1652 call errmgr%report_error(
"cholesky_factor", &
1653 "The input matrix must be square.", la_array_size_error)
1658 call dpotrf(uplo, n, a, n, flag)
1661 allocate(
character(len = 256) :: errmsg)
1662 write(errmsg, 100)
"The leading minor of order ", flag, &
1663 " is not positive definite."
1664 call errmgr%report_error(
"cholesky_factor", trim(errmsg), &
1665 la_matrix_format_error)
1669 if (uplo ==
'U')
then
1686 module subroutine cholesky_factor_cmplx(a, upper, err)
1688 complex(real64),
intent(inout),
dimension(:,:) :: a
1689 logical,
intent(in),
optional :: upper
1690 class(errors),
intent(inout),
optional,
target :: err
1693 complex(real64),
parameter :: zero = (0.0d0, 0.0d0)
1697 integer(int32) :: i, n, flag
1698 class(errors),
pointer :: errmgr
1699 type(errors),
target :: deferr
1700 character(len = :),
allocatable :: errmsg
1704 if (
present(upper))
then
1713 if (
present(err))
then
1720 if (
size(a, 2) /= n)
then
1722 call errmgr%report_error(
"cholesky_factor_cmplx", &
1723 "The input matrix must be square.", la_array_size_error)
1728 call zpotrf(uplo, n, a, n, flag)
1731 allocate(
character(len = 256) :: errmsg)
1732 write(errmsg, 100)
"The leading minor of order ", flag, &
1733 " is not positive definite."
1734 call errmgr%report_error(
"cholesky_factor_cmplx", trim(errmsg), &
1735 la_matrix_format_error)
1739 if (uplo ==
'U')
then
1756 module subroutine cholesky_rank1_update_dbl(r, u, work, err)
1758 real(real64),
intent(inout),
dimension(:,:) :: r
1759 real(real64),
intent(inout),
dimension(:) :: u
1760 real(real64),
intent(out),
target,
optional,
dimension(:) :: work
1761 class(errors),
intent(inout),
optional,
target :: err
1764 integer(int32) :: n, lwork, istat, flag
1765 real(real64),
pointer,
dimension(:) :: wptr
1766 real(real64),
allocatable,
target,
dimension(:) :: wrk
1767 class(errors),
pointer :: errmgr
1768 type(errors),
target :: deferr
1769 character(len = :),
allocatable :: errmsg
1774 if (
present(err))
then
1782 if (
size(r, 2) /= n)
then
1784 else if (
size(u) /= n)
then
1788 allocate(
character(len = 256) :: errmsg)
1789 write(errmsg, 100)
"Input number ", flag, &
1790 " is not sized correctly."
1791 call errmgr%report_error(
"cholesky_rank1_update", trim(errmsg), &
1792 la_array_size_error)
1797 if (
present(work))
then
1798 if (
size(work) < lwork)
then
1800 call errmgr%report_error(
"cholesky_rank1_update", &
1801 "The workspace array is too short.", &
1802 la_array_size_error)
1805 wptr => work(1:lwork)
1807 allocate(wrk(lwork), stat = istat)
1808 if (istat /= 0)
then
1809 call errmgr%report_error(
"cholesky_rank1_update", &
1810 "Insufficient memory available.", &
1811 la_out_of_memory_error)
1818 call dch1up(n, r, n, u, wptr)
1825 module subroutine cholesky_rank1_update_cmplx(r, u, work, err)
1827 complex(real64),
intent(inout),
dimension(:,:) :: r
1828 complex(real64),
intent(inout),
dimension(:) :: u
1829 real(real64),
intent(out),
target,
optional,
dimension(:) :: work
1830 class(errors),
intent(inout),
optional,
target :: err
1833 integer(int32) :: n, lwork, istat, flag
1834 real(real64),
pointer,
dimension(:) :: wptr
1835 real(real64),
allocatable,
target,
dimension(:) :: wrk
1836 class(errors),
pointer :: errmgr
1837 type(errors),
target :: deferr
1838 character(len = :),
allocatable :: errmsg
1843 if (
present(err))
then
1851 if (
size(r, 2) /= n)
then
1853 else if (
size(u) /= n)
then
1857 allocate(
character(len = 256) :: errmsg)
1858 write(errmsg, 100)
"Input number ", flag, &
1859 " is not sized correctly."
1860 call errmgr%report_error(
"cholesky_rank1_update_cmplx", &
1862 la_array_size_error)
1867 if (
present(work))
then
1868 if (
size(work) < lwork)
then
1870 call errmgr%report_error(
"cholesky_rank1_update_cmplx", &
1871 "The workspace array is too short.", &
1872 la_array_size_error)
1875 wptr => work(1:lwork)
1877 allocate(wrk(lwork), stat = istat)
1878 if (istat /= 0)
then
1879 call errmgr%report_error(
"cholesky_rank1_update", &
1880 "Insufficient memory available.", &
1881 la_out_of_memory_error)
1888 call zch1up(n, r, n, u, wptr)
1895 module subroutine cholesky_rank1_downdate_dbl(r, u, work, err)
1897 real(real64),
intent(inout),
dimension(:,:) :: r
1898 real(real64),
intent(inout),
dimension(:) :: u
1899 real(real64),
intent(out),
target,
optional,
dimension(:) :: work
1900 class(errors),
intent(inout),
optional,
target :: err
1903 integer(int32) :: n, lwork, istat, flag
1904 real(real64),
pointer,
dimension(:) :: wptr
1905 real(real64),
allocatable,
target,
dimension(:) :: wrk
1906 class(errors),
pointer :: errmgr
1907 type(errors),
target :: deferr
1908 character(len = :),
allocatable :: errmsg
1913 if (
present(err))
then
1921 if (
size(r, 2) /= n)
then
1923 else if (
size(u) /= n)
then
1927 allocate(
character(len = 256) :: errmsg)
1928 write(errmsg, 100)
"Input number ", flag, &
1929 " is not sized correctly."
1930 call errmgr%report_error(
"cholesky_rank1_downdate", trim(errmsg), &
1931 la_array_size_error)
1936 if (
present(work))
then
1937 if (
size(work) < lwork)
then
1939 call errmgr%report_error(
"cholesky_rank1_downdate", &
1940 "The workspace array is too short.", &
1941 la_array_size_error)
1944 wptr => work(1:lwork)
1946 allocate(wrk(lwork), stat = istat)
1947 if (istat /= 0)
then
1948 call errmgr%report_error(
"cholesky_rank1_downdate", &
1949 "Insufficient memory available.", &
1950 la_out_of_memory_error)
1957 call dch1dn(n, r, n, u, wptr, flag)
1960 call errmgr%report_error(
"cholesky_rank1_downdate", &
1961 "The downdated matrix is not positive definite.", &
1962 la_matrix_format_error)
1963 else if (flag == 2)
then
1965 call errmgr%report_error(
"cholesky_rank1_downdate", &
1966 "The input matrix is singular.", la_singular_matrix_error)
1974 module subroutine cholesky_rank1_downdate_cmplx(r, u, work, err)
1976 complex(real64),
intent(inout),
dimension(:,:) :: r
1977 complex(real64),
intent(inout),
dimension(:) :: u
1978 real(real64),
intent(out),
target,
optional,
dimension(:) :: work
1979 class(errors),
intent(inout),
optional,
target :: err
1982 integer(int32) :: n, lwork, istat, flag
1983 real(real64),
pointer,
dimension(:) :: wptr
1984 real(real64),
allocatable,
target,
dimension(:) :: wrk
1985 class(errors),
pointer :: errmgr
1986 type(errors),
target :: deferr
1987 character(len = :),
allocatable :: errmsg
1992 if (
present(err))
then
2000 if (
size(r, 2) /= n)
then
2002 else if (
size(u) /= n)
then
2006 allocate(
character(len = 256) :: errmsg)
2007 write(errmsg, 100)
"Input number ", flag, &
2008 " is not sized correctly."
2009 call errmgr%report_error(
"cholesky_rank1_downdate_cmplx", &
2011 la_array_size_error)
2016 if (
present(work))
then
2017 if (
size(work) < lwork)
then
2019 call errmgr%report_error(
"cholesky_rank1_downdate_cmplx", &
2020 "The workspace array is too short.", &
2021 la_array_size_error)
2024 wptr => work(1:lwork)
2026 allocate(wrk(lwork), stat = istat)
2027 if (istat /= 0)
then
2028 call errmgr%report_error(
"cholesky_rank1_downdate_cmplx", &
2029 "Insufficient memory available.", &
2030 la_out_of_memory_error)
2037 call zch1dn(n, r, n, u, wptr, flag)
2040 call errmgr%report_error(
"cholesky_rank1_downdate_cmplx", &
2041 "The downdated matrix is not positive definite.", &
2042 la_matrix_format_error)
2043 else if (flag == 2)
then
2045 call errmgr%report_error(
"cholesky_rank1_downdate_cmplx", &
2046 "The input matrix is singular.", la_singular_matrix_error)
2056 module subroutine rz_factor_dbl(a, tau, work, olwork, err)
2058 real(real64),
intent(inout),
dimension(:,:) :: a
2059 real(real64),
intent(out),
dimension(:) :: tau
2060 real(real64),
intent(out),
target,
optional,
dimension(:) :: work
2061 integer(int32),
intent(out),
optional :: olwork
2062 class(errors),
intent(inout),
optional,
target :: err
2065 integer(int32) :: m, n, lwork, flag, istat
2066 real(real64),
pointer,
dimension(:) :: wptr
2067 real(real64),
allocatable,
target,
dimension(:) :: wrk
2068 real(real64),
dimension(1) :: temp
2069 class(errors),
pointer :: errmgr
2070 type(errors),
target :: deferr
2071 character(len = :),
allocatable :: errmsg
2076 if (
present(err))
then
2084 if (
size(tau) /= m)
then
2089 allocate(
character(len = 256) :: errmsg)
2090 write(errmsg, 100)
"Input number ", flag, &
2091 " is not sized correctly."
2092 call errmgr%report_error(
"rz_factor", trim(errmsg), &
2093 la_array_size_error)
2098 call dtzrzf(m, n, a, m, tau, temp, -1, flag)
2099 lwork = int(temp(1), int32)
2100 if (
present(olwork))
then
2106 if (
present(work))
then
2107 if (
size(work) < lwork)
then
2109 call errmgr%report_error(
"rz_factor", &
2110 "Incorrectly sized input array WORK, argument 3.", &
2111 la_array_size_error)
2114 wptr => work(1:lwork)
2116 allocate(wrk(lwork), stat = istat)
2117 if (istat /= 0)
then
2119 call errmgr%report_error(
"rz_factor", &
2120 "Insufficient memory available.", &
2121 la_out_of_memory_error)
2128 call dtzrzf(m, n, a, m, tau, wptr, lwork, flag)
2135 module subroutine rz_factor_cmplx(a, tau, work, olwork, err)
2137 complex(real64),
intent(inout),
dimension(:,:) :: a
2138 complex(real64),
intent(out),
dimension(:) :: tau
2139 complex(real64),
intent(out),
target,
optional,
dimension(:) :: work
2140 integer(int32),
intent(out),
optional :: olwork
2141 class(errors),
intent(inout),
optional,
target :: err
2144 integer(int32) :: m, n, lwork, flag, istat
2145 complex(real64),
pointer,
dimension(:) :: wptr
2146 complex(real64),
allocatable,
target,
dimension(:) :: wrk
2147 complex(real64),
dimension(1) :: temp
2148 class(errors),
pointer :: errmgr
2149 type(errors),
target :: deferr
2150 character(len = :),
allocatable :: errmsg
2155 if (
present(err))
then
2163 if (
size(tau) /= m)
then
2168 allocate(
character(len = 256) :: errmsg)
2169 write(errmsg, 100)
"Input number ", flag, &
2170 " is not sized correctly."
2171 call errmgr%report_error(
"rz_factor_cmplx", trim(errmsg), &
2172 la_array_size_error)
2177 call ztzrzf(m, n, a, m, tau, temp, -1, flag)
2178 lwork = int(temp(1), int32)
2179 if (
present(olwork))
then
2185 if (
present(work))
then
2186 if (
size(work) < lwork)
then
2188 call errmgr%report_error(
"rz_factor_cmplx", &
2189 "Incorrectly sized input array WORK, argument 3.", &
2190 la_array_size_error)
2193 wptr => work(1:lwork)
2195 allocate(wrk(lwork), stat = istat)
2196 if (istat /= 0)
then
2198 call errmgr%report_error(
"rz_factor_cmplx", &
2199 "Insufficient memory available.", &
2200 la_out_of_memory_error)
2207 call ztzrzf(m, n, a, m, tau, wptr, lwork, flag)
2214 module subroutine mult_rz_mtx(lside, trans, l, a, tau, c, work, olwork, err)
2216 logical,
intent(in) :: lside, trans
2217 integer(int32),
intent(in) :: l
2218 real(real64),
intent(inout),
dimension(:,:) :: a, c
2219 real(real64),
intent(in),
dimension(:) :: tau
2220 real(real64),
intent(out),
target,
optional,
dimension(:) :: work
2221 integer(int32),
intent(out),
optional :: olwork
2222 class(errors),
intent(inout),
optional,
target :: err
2225 character :: side, t
2226 integer(int32) :: m, n, k, lwork, flag, istat, lda
2227 real(real64),
pointer,
dimension(:) :: wptr
2228 real(real64),
allocatable,
target,
dimension(:) :: wrk
2229 real(real64),
dimension(1) :: temp
2230 class(errors),
pointer :: errmgr
2231 type(errors),
target :: deferr
2232 character(len = :),
allocatable :: errmsg
2249 if (
present(err))
then
2258 if (l > m .or. l < 0)
then
2260 else if (k > m)
then
2262 else if (
size(a, 1) < k .or.
size(a, 2) /= m)
then
2266 if (l > n .or. l < 0)
then
2268 else if (k > n)
then
2270 else if (
size(a, 1) < k .or.
size(a, 2) /= n)
then
2276 allocate(
character(len = 256) :: errmsg)
2277 write(errmsg, 100)
"Input number ", flag, &
2278 " is not sized correctly."
2279 call errmgr%report_error(
"mult_rz_mtx", trim(errmsg), &
2280 la_array_size_error)
2285 call dormrz(side, t, m, n, k, l, a, lda, tau, c, m, temp, -1, flag)
2286 lwork = int(temp(1), int32)
2287 if (
present(olwork))
then
2293 if (
present(work))
then
2294 if (
size(work) < lwork)
then
2296 call errmgr%report_error(
"mult_rz_mtx", &
2297 "Incorrectly sized input array WORK, argument 7.", &
2298 la_array_size_error)
2301 wptr => work(1:lwork)
2303 allocate(wrk(lwork), stat = istat)
2304 if (istat /= 0)
then
2306 call errmgr%report_error(
"mult_rz_mtx", &
2307 "Insufficient memory available.", &
2308 la_out_of_memory_error)
2315 call dormrz(side, t, m, n, k, l, a, lda, tau, c, m, wptr, lwork, flag)
2322 module subroutine mult_rz_mtx_cmplx(lside, trans, l, a, tau, c, work, olwork, err)
2324 logical,
intent(in) :: lside, trans
2325 integer(int32),
intent(in) :: l
2326 complex(real64),
intent(inout),
dimension(:,:) :: a, c
2327 complex(real64),
intent(in),
dimension(:) :: tau
2328 complex(real64),
intent(out),
target,
optional,
dimension(:) :: work
2329 integer(int32),
intent(out),
optional :: olwork
2330 class(errors),
intent(inout),
optional,
target :: err
2333 character :: side, t
2334 integer(int32) :: m, n, k, lwork, flag, istat, lda
2335 complex(real64),
pointer,
dimension(:) :: wptr
2336 complex(real64),
allocatable,
target,
dimension(:) :: wrk
2337 complex(real64),
dimension(1) :: temp
2338 class(errors),
pointer :: errmgr
2339 type(errors),
target :: deferr
2340 character(len = :),
allocatable :: errmsg
2357 if (
present(err))
then
2366 if (l > m .or. l < 0)
then
2368 else if (k > m)
then
2370 else if (
size(a, 1) < k .or.
size(a, 2) /= m)
then
2374 if (l > n .or. l < 0)
then
2376 else if (k > n)
then
2378 else if (
size(a, 1) < k .or.
size(a, 2) /= n)
then
2384 allocate(
character(len = 256) :: errmsg)
2385 write(errmsg, 100)
"Input number ", flag, &
2386 " is not sized correctly."
2387 call errmgr%report_error(
"mult_rz_mtx_cmplx", trim(errmsg), &
2388 la_array_size_error)
2393 call zunmrz(side, t, m, n, k, l, a, lda, tau, c, m, temp, -1, flag)
2394 lwork = int(temp(1), int32)
2395 if (
present(olwork))
then
2401 if (
present(work))
then
2402 if (
size(work) < lwork)
then
2404 call errmgr%report_error(
"mult_rz_mtx_cmplx", &
2405 "Incorrectly sized input array WORK, argument 7.", &
2406 la_array_size_error)
2409 wptr => work(1:lwork)
2411 allocate(wrk(lwork), stat = istat)
2412 if (istat /= 0)
then
2414 call errmgr%report_error(
"mult_rz_mtx_cmplx", &
2415 "Insufficient memory available.", &
2416 la_out_of_memory_error)
2423 call zunmrz(side, t, m, n, k, l, a, lda, tau, c, m, wptr, lwork, flag)
2430 module subroutine mult_rz_vec(trans, l, a, tau, c, work, olwork, err)
2432 logical,
intent(in) :: trans
2433 integer(int32),
intent(in) :: l
2434 real(real64),
intent(inout),
dimension(:,:) :: a
2435 real(real64),
intent(in),
dimension(:) :: tau
2436 real(real64),
intent(inout),
dimension(:) :: c
2437 real(real64),
intent(out),
target,
optional,
dimension(:) :: work
2438 integer(int32),
intent(out),
optional :: olwork
2439 class(errors),
intent(inout),
optional,
target :: err
2442 character :: side, t
2443 integer(int32) :: m, k, lwork, flag, istat, lda
2444 real(real64),
pointer,
dimension(:) :: wptr
2445 real(real64),
allocatable,
target,
dimension(:) :: wrk
2446 real(real64),
dimension(1) :: temp
2447 class(errors),
pointer :: errmgr
2448 type(errors),
target :: deferr
2449 character(len = :),
allocatable :: errmsg
2461 if (
present(err))
then
2469 if (l > m .or. l < 0)
then
2471 else if (k > m)
then
2473 else if (
size(a, 1) < k .or.
size(a, 2) /= m)
then
2478 allocate(
character(len = 256) :: errmsg)
2479 write(errmsg, 100)
"Input number ", flag, &
2480 " is not sized correctly."
2481 call errmgr%report_error(
"mult_rz_vec", trim(errmsg), &
2482 la_array_size_error)
2487 call dormrz(side, t, m, 1, k, l, a, lda, tau, c, m, temp, -1, flag)
2488 lwork = int(temp(1), int32)
2489 if (
present(olwork))
then
2495 if (
present(work))
then
2496 if (
size(work) < lwork)
then
2498 call errmgr%report_error(
"mult_rz_vec", &
2499 "Incorrectly sized input array WORK, argument 6.", &
2500 la_array_size_error)
2503 wptr => work(1:lwork)
2505 allocate(wrk(lwork), stat = istat)
2506 if (istat /= 0)
then
2508 call errmgr%report_error(
"mult_rz_vec", &
2509 "Insufficient memory available.", &
2510 la_out_of_memory_error)
2517 call dormrz(side, t, m, 1, k, l, a, lda, tau, c, m, wptr, lwork, flag)
2524 module subroutine mult_rz_vec_cmplx(trans, l, a, tau, c, work, olwork, err)
2526 logical,
intent(in) :: trans
2527 integer(int32),
intent(in) :: l
2528 complex(real64),
intent(inout),
dimension(:,:) :: a
2529 complex(real64),
intent(in),
dimension(:) :: tau
2530 complex(real64),
intent(inout),
dimension(:) :: c
2531 complex(real64),
intent(out),
target,
optional,
dimension(:) :: work
2532 integer(int32),
intent(out),
optional :: olwork
2533 class(errors),
intent(inout),
optional,
target :: err
2536 character :: side, t
2537 integer(int32) :: m, k, lwork, flag, istat, lda
2538 complex(real64),
pointer,
dimension(:) :: wptr
2539 complex(real64),
allocatable,
target,
dimension(:) :: wrk
2540 complex(real64),
dimension(1) :: temp
2541 class(errors),
pointer :: errmgr
2542 type(errors),
target :: deferr
2543 character(len = :),
allocatable :: errmsg
2555 if (
present(err))
then
2563 if (l > m .or. l < 0)
then
2565 else if (k > m)
then
2567 else if (
size(a, 1) < k .or.
size(a, 2) /= m)
then
2572 allocate(
character(len = 256) :: errmsg)
2573 write(errmsg, 100)
"Input number ", flag, &
2574 " is not sized correctly."
2575 call errmgr%report_error(
"mult_rz_vec_cmplx", trim(errmsg), &
2576 la_array_size_error)
2581 call zunmrz(side, t, m, 1, k, l, a, lda, tau, c, m, temp, -1, flag)
2582 lwork = int(temp(1), int32)
2583 if (
present(olwork))
then
2589 if (
present(work))
then
2590 if (
size(work) < lwork)
then
2592 call errmgr%report_error(
"mult_rz_vec_cmplx", &
2593 "Incorrectly sized input array WORK, argument 6.", &
2594 la_array_size_error)
2597 wptr => work(1:lwork)
2599 allocate(wrk(lwork), stat = istat)
2600 if (istat /= 0)
then
2602 call errmgr%report_error(
"mult_rz_vec_cmplx", &
2603 "Insufficient memory available.", &
2604 la_out_of_memory_error)
2611 call zunmrz(side, t, m, 1, k, l, a, lda, tau, c, m, wptr, lwork, flag)
2620 module subroutine svd_dbl(a, s, u, vt, work, olwork, err)
2622 real(real64),
intent(inout),
dimension(:,:) :: a
2623 real(real64),
intent(out),
dimension(:) :: s
2624 real(real64),
intent(out),
optional,
dimension(:,:) :: u, vt
2625 real(real64),
intent(out),
target,
optional,
dimension(:) :: work
2626 integer(int32),
intent(out),
optional :: olwork
2627 class(errors),
intent(inout),
optional,
target :: err
2630 character :: jobu, jobvt
2631 integer(int32) :: m, n, mn, istat, lwork, flag
2632 real(real64),
pointer,
dimension(:) :: wptr
2633 real(real64),
allocatable,
target,
dimension(:) :: wrk
2634 real(real64),
dimension(1) :: temp
2635 class(errors),
pointer :: errmgr
2636 type(errors),
target :: deferr
2637 character(len = :),
allocatable :: errmsg
2643 if (
present(u))
then
2644 if (
size(u, 2) == m)
then
2646 else if (
size(u, 2) == mn)
then
2652 if (
present(vt))
then
2657 if (
present(err))
then
2665 if (
size(s) /= mn)
then
2667 else if (
present(u))
then
2668 if (
size(u, 1) /= m) flag = 3
2669 if (
size(u, 2) /= m .and.
size(u, 2) /= mn) flag = 3
2670 else if (
present(vt))
then
2671 if (
size(vt, 1) /= n .or.
size(vt, 2) /= n) flag = 4
2675 allocate(
character(len = 256) :: errmsg)
2676 write(errmsg, 100)
"Input number ", flag, &
2677 " is not sized correctly."
2678 call errmgr%report_error(
"svd", trim(errmsg), &
2679 la_array_size_error)
2684 call dgesvd(jobu, jobvt, m, n, a, m, s, temp, m, temp, n, temp, -1, &
2686 lwork = int(temp(1), int32)
2687 if (
present(olwork))
then
2693 if (
present(work))
then
2694 if (
size(work) < lwork)
then
2696 call errmgr%report_error(
"svd", &
2697 "Incorrectly sized input array WORK, argument 5.", &
2698 la_array_size_error)
2701 wptr => work(1:lwork)
2703 allocate(wrk(lwork), stat = istat)
2704 if (istat /= 0)
then
2706 call errmgr%report_error(
"svd", &
2707 "Insufficient memory available.", &
2708 la_out_of_memory_error)
2715 if (
present(u) .and.
present(vt))
then
2716 call dgesvd(jobu, jobvt, m, n, a, m, s, u, m, vt, n, wptr, lwork, &
2718 else if (
present(u) .and. .not.
present(vt))
then
2719 call dgesvd(jobu, jobvt, m, n, a, m, s, u, m, temp, n, wptr, &
2721 else if (.not.
present(u) .and.
present(vt))
then
2722 call dgesvd(jobu, jobvt, m, n, a, m, s, temp, m, vt, n, wptr, &
2725 call dgesvd(jobu, jobvt, m, n, a, m, s, temp, m, temp, n, wptr, &
2731 allocate(
character(len = 256) :: errmsg)
2732 write(errmsg, 101) flag,
" superdiagonals could not " // &
2733 "converge to zero as part of the QR iteration process."
2734 call errmgr%report_warning(
"svd", errmsg, la_convergence_error)
2743 module subroutine svd_cmplx(a, s, u, vt, work, olwork, rwork, err)
2745 complex(real64),
intent(inout),
dimension(:,:) :: a
2746 real(real64),
intent(out),
dimension(:) :: s
2747 complex(real64),
intent(out),
optional,
dimension(:,:) :: u, vt
2748 complex(real64),
intent(out),
target,
optional,
dimension(:) :: work
2749 integer(int32),
intent(out),
optional :: olwork
2750 real(real64),
intent(out),
target,
optional,
dimension(:) :: rwork
2751 class(errors),
intent(inout),
optional,
target :: err
2754 character :: jobu, jobvt
2755 integer(int32) :: m, n, mn, istat, lwork, flag, lrwork
2756 complex(real64),
pointer,
dimension(:) :: wptr
2757 complex(real64),
allocatable,
target,
dimension(:) :: wrk
2758 complex(real64),
dimension(1) :: temp
2759 real(real64),
dimension(1) :: rtemp
2760 real(real64),
pointer,
dimension(:) :: rwptr
2761 real(real64),
allocatable,
target,
dimension(:) :: rwrk
2762 class(errors),
pointer :: errmgr
2763 type(errors),
target :: deferr
2764 character(len = :),
allocatable :: errmsg
2771 if (
present(u))
then
2772 if (
size(u, 2) == m)
then
2774 else if (
size(u, 2) == mn)
then
2780 if (
present(vt))
then
2785 if (
present(err))
then
2793 if (
size(s) /= mn)
then
2795 else if (
present(u))
then
2796 if (
size(u, 1) /= m) flag = 3
2797 if (
size(u, 2) /= m .and.
size(u, 2) /= mn) flag = 3
2798 else if (
present(vt))
then
2799 if (
size(vt, 1) /= n .or.
size(vt, 2) /= n) flag = 4
2803 allocate(
character(len = 256) :: errmsg)
2804 write(errmsg, 100)
"Input number ", flag, &
2805 " is not sized correctly."
2806 call errmgr%report_error(
"svd_cmplx", trim(errmsg), &
2807 la_array_size_error)
2812 call zgesvd(jobu, jobvt, m, n, a, m, s, temp, m, temp, n, temp, -1, &
2814 lwork = int(temp(1), int32)
2815 if (
present(olwork))
then
2821 if (
present(work))
then
2822 if (
size(work) < lwork)
then
2824 call errmgr%report_error(
"svd_cmplx", &
2825 "Incorrectly sized input array WORK, argument 5.", &
2826 la_array_size_error)
2829 wptr => work(1:lwork)
2831 allocate(wrk(lwork), stat = istat)
2832 if (istat /= 0)
then
2834 call errmgr%report_error(
"svd_cmplx", &
2835 "Insufficient memory available.", &
2836 la_out_of_memory_error)
2842 if (
present(rwork))
then
2843 if (
size(rwork) < lrwork)
then
2845 call errmgr%report_error(
"svd_cmplx", &
2846 "Incorrectly sized input array RWORK, argument 7.", &
2847 la_array_size_error)
2849 rwptr => rwork(1:lrwork)
2851 allocate(rwrk(lrwork), stat = istat)
2852 if (istat /= 0)
then
2854 call errmgr%report_error(
"svd_cmplx", &
2855 "Insufficient memory available.", &
2856 la_out_of_memory_error)
2863 if (
present(u) .and.
present(vt))
then
2864 call zgesvd(jobu, jobvt, m, n, a, m, s, u, m, vt, n, wptr, lwork, &
2866 else if (
present(u) .and. .not.
present(vt))
then
2867 call zgesvd(jobu, jobvt, m, n, a, m, s, u, m, temp, n, wptr, &
2869 else if (.not.
present(u) .and.
present(vt))
then
2870 call zgesvd(jobu, jobvt, m, n, a, m, s, temp, m, vt, n, wptr, &
2873 call zgesvd(jobu, jobvt, m, n, a, m, s, temp, m, temp, n, wptr, &
2879 allocate(
character(len = 256) :: errmsg)
2880 write(errmsg, 101) flag,
" superdiagonals could not " // &
2881 "converge to zero as part of the QR iteration process."
2882 call errmgr%report_warning(
"svd_cmplx", errmsg, &
2883 la_convergence_error)
2894 module subroutine lq_factor_no_pivot(a, tau, work, olwork, err)
2896 real(real64),
intent(inout),
dimension(:,:) :: a
2897 real(real64),
intent(out),
dimension(:) :: tau
2898 real(real64),
intent(out),
target,
dimension(:),
optional :: work
2899 integer(int32),
intent(out),
optional :: olwork
2900 class(errors),
intent(inout),
optional,
target :: err
2903 integer(int32) :: m, n, mn, istat, lwork, flag
2904 real(real64),
dimension(1) :: temp
2905 real(real64),
pointer,
dimension(:) :: wptr
2906 real(real64),
allocatable,
target,
dimension(:) :: wrk
2907 class(errors),
pointer :: errmgr
2908 type(errors),
target :: deferr
2914 if (
present(err))
then
2921 if (
size(tau) /= mn)
then
2923 call errmgr%report_error(
"lq_factor_no_pivot", &
2924 "Incorrectly sized input array TAU, argument 2.", &
2925 la_array_size_error)
2930 call dgelqf(m, n, a, m, tau, temp, -1, flag)
2931 lwork = int(temp(1), int32)
2932 if (
present(olwork))
then
2938 if (
present(work))
then
2939 if (
size(work) < lwork)
then
2941 call errmgr%report_error(
"lq_factor_no_pivot", &
2942 "Incorrectly sized input array WORK, argument 3.", &
2943 la_array_size_error)
2946 wptr => work(1:lwork)
2948 allocate(wrk(lwork), stat = istat)
2949 if (istat /= 0)
then
2951 call errmgr%report_error(
"lq_factor_no_pivot", &
2952 "Insufficient memory available.", &
2953 la_out_of_memory_error)
2960 call dgelqf(m, n, a, m, tau, wptr, lwork, flag)
2964 module subroutine lq_factor_no_pivot_cmplx(a, tau, work, olwork, err)
2966 complex(real64),
intent(inout),
dimension(:,:) :: a
2967 complex(real64),
intent(out),
dimension(:) :: tau
2968 complex(real64),
intent(out),
target,
dimension(:),
optional :: work
2969 integer(int32),
intent(out),
optional :: olwork
2970 class(errors),
intent(inout),
optional,
target :: err
2973 integer(int32) :: m, n, mn, istat, lwork, flag
2974 complex(real64),
dimension(1) :: temp
2975 complex(real64),
pointer,
dimension(:) :: wptr
2976 complex(real64),
allocatable,
target,
dimension(:) :: wrk
2977 class(errors),
pointer :: errmgr
2978 type(errors),
target :: deferr
2984 if (
present(err))
then
2991 if (
size(tau) /= mn)
then
2993 call errmgr%report_error(
"lq_factor_no_pivot_cmplx", &
2994 "Incorrectly sized input array TAU, argument 2.", &
2995 la_array_size_error)
3000 call zgelqf(m, n, a, m, tau, temp, -1, flag)
3001 lwork = int(temp(1), int32)
3002 if (
present(olwork))
then
3008 if (
present(work))
then
3009 if (
size(work) < lwork)
then
3011 call errmgr%report_error(
"lq_factor_no_pivot_cmplx", &
3012 "Incorrectly sized input array WORK, argument 3.", &
3013 la_array_size_error)
3016 wptr => work(1:lwork)
3018 allocate(wrk(lwork), stat = istat)
3019 if (istat /= 0)
then
3021 call errmgr%report_error(
"lq_factor_no_pivot_cmplx", &
3022 "Insufficient memory available.", &
3023 la_out_of_memory_error)
3030 call zgelqf(m, n, a, m, tau, wptr, lwork, flag)
3034 module subroutine form_lq_no_pivot(l, tau, q, work, olwork, err)
3036 real(real64),
intent(inout),
dimension(:,:) :: l
3037 real(real64),
intent(in),
dimension(:) :: tau
3038 real(real64),
intent(out),
dimension(:,:) :: q
3039 real(real64),
intent(out),
target,
dimension(:),
optional :: work
3040 integer(int32),
intent(out),
optional :: olwork
3041 class(errors),
intent(inout),
optional,
target :: err
3044 real(real64),
parameter :: zero = 0.0d0
3047 integer(int32) :: i, m, n, mn, k, istat, flag, lwork
3048 real(real64),
pointer,
dimension(:) :: wptr
3049 real(real64),
allocatable,
target,
dimension(:) :: wrk
3050 real(real64),
dimension(1) :: temp
3051 class(errors),
pointer :: errmgr
3052 type(errors),
target :: deferr
3053 character(len = :),
allocatable :: errmsg
3059 if (
present(err))
then
3069 else if (
size(tau) /= mn)
then
3071 else if (
size(q, 1) /= n .or.
size(q, 2) /= n)
then
3076 allocate(
character(len = 256) :: errmsg)
3077 write(errmsg, 100)
"Input number ", flag, &
3078 " is not sized correctly."
3079 call errmgr%report_error(
"form_lq_no_pivot", trim(errmsg), &
3080 la_array_size_error)
3085 call dorglq(n, n, mn, q, n, tau, temp, -1, flag)
3086 lwork = int(temp(1), int32)
3087 if (
present(olwork))
then
3093 if (
present(work))
then
3094 if (
size(work) < lwork)
then
3096 call errmgr%report_error(
"form_lq_no_pivot", &
3097 "Incorrectly sized input array WORK, argument 4.", &
3098 la_array_size_error)
3101 wptr => work(1:lwork)
3103 allocate(wrk(lwork), stat = istat)
3104 if (istat /= 0)
then
3106 call errmgr%report_error(
"form_lq_no_pivot", &
3107 "Insufficient memory available.", &
3108 la_out_of_memory_error)
3117 q(1:j-1,j) = l(1:k,j)
3122 call dorglq(n, n, mn, q, n, tau, wptr, lwork, flag)
3129 module subroutine form_lq_no_pivot_cmplx(l, tau, q, work, olwork, err)
3131 complex(real64),
intent(inout),
dimension(:,:) :: l
3132 complex(real64),
intent(in),
dimension(:) :: tau
3133 complex(real64),
intent(out),
dimension(:,:) :: q
3134 complex(real64),
intent(out),
target,
dimension(:),
optional :: work
3135 integer(int32),
intent(out),
optional :: olwork
3136 class(errors),
intent(inout),
optional,
target :: err
3139 complex(real64),
parameter :: zero = (0.0d0, 0.0d0)
3142 integer(int32) :: i, m, n, mn, k, istat, flag, lwork
3143 complex(real64),
pointer,
dimension(:) :: wptr
3144 complex(real64),
allocatable,
target,
dimension(:) :: wrk
3145 complex(real64),
dimension(1) :: temp
3146 class(errors),
pointer :: errmgr
3147 type(errors),
target :: deferr
3148 character(len = :),
allocatable :: errmsg
3154 if (
present(err))
then
3164 else if (
size(tau) /= mn)
then
3166 else if (
size(q, 1) /= n .or.
size(q, 2) /= n)
then
3171 allocate(
character(len = 256) :: errmsg)
3172 write(errmsg, 100)
"Input number ", flag, &
3173 " is not sized correctly."
3174 call errmgr%report_error(
"form_lq_no_pivot_cmplx", trim(errmsg), &
3175 la_array_size_error)
3180 call zunglq(n, n, mn, q, n, tau, temp, -1, flag)
3181 lwork = int(temp(1), int32)
3182 if (
present(olwork))
then
3188 if (
present(work))
then
3189 if (
size(work) < lwork)
then
3191 call errmgr%report_error(
"form_lq_no_pivot_cmplx", &
3192 "Incorrectly sized input array WORK, argument 4.", &
3193 la_array_size_error)
3196 wptr => work(1:lwork)
3198 allocate(wrk(lwork), stat = istat)
3199 if (istat /= 0)
then
3201 call errmgr%report_error(
"form_lq_no_pivot_cmplx", &
3202 "Insufficient memory available.", &
3203 la_out_of_memory_error)
3212 q(1:j-1,j) = l(1:k,j)
3217 call zunglq(n, n, mn, q, n, tau, wptr, lwork, flag)
3224 module subroutine mult_lq_mtx(lside, trans, a, tau, c, work, olwork, err)
3226 logical,
intent(in) :: lside, trans
3227 real(real64),
intent(in),
dimension(:,:) :: a
3228 real(real64),
intent(in),
dimension(:) :: tau
3229 real(real64),
intent(inout),
dimension(:,:) :: c
3230 real(real64),
intent(out),
target,
dimension(:),
optional :: work
3231 integer(int32),
intent(out),
optional :: olwork
3232 class(errors),
intent(inout),
optional,
target :: err
3235 character :: side, t
3236 integer(int32) :: m, n, k, ncola, istat, flag, lwork
3237 real(real64),
pointer,
dimension(:) :: wptr
3238 real(real64),
allocatable,
target,
dimension(:) :: wrk
3239 real(real64),
dimension(1) :: temp
3240 class(errors),
pointer :: errmgr
3241 type(errors),
target :: deferr
3242 character(len = :),
allocatable :: errmsg
3260 if (
present(err))
then
3268 if (
size(a, 1) /= k .or.
size(a, 2) /= ncola)
then
3273 allocate(
character(len = 256) :: errmsg)
3274 write(errmsg, 100)
"Input number ", flag, &
3275 " is not sized correctly."
3276 call errmgr%report_error(
"mult_lq_mtx", trim(errmsg), &
3277 la_array_size_error)
3282 call dormlq(side, t, m, n, k, a, k, tau, c, m, temp, -1, flag)
3283 lwork = int(temp(1), int32)
3284 if (
present(olwork))
then
3290 if (
present(work))
then
3291 if (
size(work) < lwork)
then
3293 call errmgr%report_error(
"mult_lq_mtx", &
3294 "Incorrectly sized input array WORK, argument 6.", &
3295 la_array_size_error)
3298 wptr => work(1:lwork)
3300 allocate(wrk(lwork), stat = istat)
3301 if (istat /= 0)
then
3303 call errmgr%report_error(
"mult_lq_mtx", &
3304 "Insufficient memory available.", &
3305 la_out_of_memory_error)
3312 call dormlq(side, t, m, n, k, a, k, tau, c, m, wptr, lwork, flag)
3319 module subroutine mult_lq_mtx_cmplx(lside, trans, a, tau, c, work, olwork, err)
3321 logical,
intent(in) :: lside, trans
3322 complex(real64),
intent(in),
dimension(:,:) :: a
3323 complex(real64),
intent(in),
dimension(:) :: tau
3324 complex(real64),
intent(inout),
dimension(:,:) :: c
3325 complex(real64),
intent(out),
target,
dimension(:),
optional :: work
3326 integer(int32),
intent(out),
optional :: olwork
3327 class(errors),
intent(inout),
optional,
target :: err
3330 character :: side, t
3331 integer(int32) :: m, n, k, ncola, istat, flag, lwork
3332 complex(real64),
pointer,
dimension(:) :: wptr
3333 complex(real64),
allocatable,
target,
dimension(:) :: wrk
3334 complex(real64),
dimension(1) :: temp
3335 class(errors),
pointer :: errmgr
3336 type(errors),
target :: deferr
3337 character(len = :),
allocatable :: errmsg
3355 if (
present(err))
then
3363 if (
size(a, 1) /= k .or.
size(a, 2) /= ncola)
then
3368 allocate(
character(len = 256) :: errmsg)
3369 write(errmsg, 100)
"Input number ", flag, &
3370 " is not sized correctly."
3371 call errmgr%report_error(
"mult_lq_mtx_cmplx", trim(errmsg), &
3372 la_array_size_error)
3377 call zunmlq(side, t, m, n, k, a, k, tau, c, m, temp, -1, flag)
3378 lwork = int(temp(1), int32)
3379 if (
present(olwork))
then
3385 if (
present(work))
then
3386 if (
size(work) < lwork)
then
3388 call errmgr%report_error(
"mult_lq_mtx_cmplx", &
3389 "Incorrectly sized input array WORK, argument 6.", &
3390 la_array_size_error)
3393 wptr => work(1:lwork)
3395 allocate(wrk(lwork), stat = istat)
3396 if (istat /= 0)
then
3398 call errmgr%report_error(
"mult_lq_mtx_cmplx", &
3399 "Insufficient memory available.", &
3400 la_out_of_memory_error)
3407 call zunmlq(side, t, m, n, k, a, k, tau, c, m, wptr, lwork, flag)
3414 module subroutine mult_lq_vec(trans, a, tau, c, work, olwork, err)
3416 logical,
intent(in) :: trans
3417 real(real64),
intent(in),
dimension(:,:) :: a
3418 real(real64),
intent(in),
dimension(:) :: tau
3419 real(real64),
intent(inout),
dimension(:) :: c
3420 real(real64),
intent(out),
target,
dimension(:),
optional :: work
3421 integer(int32),
intent(out),
optional :: olwork
3422 class(errors),
intent(inout),
optional,
target :: err
3425 character :: side, t
3426 integer(int32) :: m, n, k, istat, flag, lwork
3427 real(real64),
pointer,
dimension(:) :: wptr
3428 real(real64),
allocatable,
target,
dimension(:) :: wrk
3429 real(real64),
dimension(1) :: temp
3430 class(errors),
pointer :: errmgr
3431 type(errors),
target :: deferr
3432 character(len = :),
allocatable :: errmsg
3444 if (
present(err))
then
3452 if (
size(a, 1) /= k .or.
size(a, 2) /= m)
then
3457 allocate(
character(len = 256) :: errmsg)
3458 write(errmsg, 100)
"Input number ", flag, &
3459 " is not sized correctly."
3460 call errmgr%report_error(
"mult_lq_vec", trim(errmsg), &
3461 la_array_size_error)
3466 call dormlq(side, t, m, n, k, a, k, tau, c, m, temp, -1, flag)
3467 lwork = int(temp(1), int32)
3468 if (
present(olwork))
then
3474 if (
present(work))
then
3475 if (
size(work) < lwork)
then
3477 call errmgr%report_error(
"mult_lq_vec", &
3478 "Incorrectly sized input array WORK, argument 6.", &
3479 la_array_size_error)
3482 wptr => work(1:lwork)
3484 allocate(wrk(lwork), stat = istat)
3485 if (istat /= 0)
then
3487 call errmgr%report_error(
"mult_lq_vec", &
3488 "Insufficient memory available.", &
3489 la_out_of_memory_error)
3496 call dormlq(side, t, m, n, k, a, k, tau, c, m, wptr, lwork, flag)
3503 module subroutine mult_lq_vec_cmplx(trans, a, tau, c, work, olwork, err)
3505 logical,
intent(in) :: trans
3506 complex(real64),
intent(in),
dimension(:,:) :: a
3507 complex(real64),
intent(in),
dimension(:) :: tau
3508 complex(real64),
intent(inout),
dimension(:) :: c
3509 complex(real64),
intent(out),
target,
dimension(:),
optional :: work
3510 integer(int32),
intent(out),
optional :: olwork
3511 class(errors),
intent(inout),
optional,
target :: err
3514 character :: side, t
3515 integer(int32) :: m, n, k, istat, flag, lwork
3516 complex(real64),
pointer,
dimension(:) :: wptr
3517 complex(real64),
allocatable,
target,
dimension(:) :: wrk
3518 complex(real64),
dimension(1) :: temp
3519 class(errors),
pointer :: errmgr
3520 type(errors),
target :: deferr
3521 character(len = :),
allocatable :: errmsg
3533 if (
present(err))
then
3541 if (
size(a, 1) /= k .or.
size(a, 2) /= m)
then
3546 allocate(
character(len = 256) :: errmsg)
3547 write(errmsg, 100)
"Input number ", flag, &
3548 " is not sized correctly."
3549 call errmgr%report_error(
"mult_lq_vec_cmplx", trim(errmsg), &
3550 la_array_size_error)
3555 call zunmlq(side, t, m, n, k, a, k, tau, c, m, temp, -1, flag)
3556 lwork = int(temp(1), int32)
3557 if (
present(olwork))
then
3563 if (
present(work))
then
3564 if (
size(work) < lwork)
then
3566 call errmgr%report_error(
"mult_lq_vec_cmplx", &
3567 "Incorrectly sized input array WORK, argument 6.", &
3568 la_array_size_error)
3571 wptr => work(1:lwork)
3573 allocate(wrk(lwork), stat = istat)
3574 if (istat /= 0)
then
3576 call errmgr%report_error(
"mult_lq_vec_cmplx", &
3577 "Insufficient memory available.", &
3578 la_out_of_memory_error)
3585 call zunmlq(side, t, m, n, k, a, k, tau, c, m, wptr, lwork, flag)
Provides a set of common linear algebra routines.